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FOREWORD 

The  analyses  described  in  this  report  were  performed  by  TRW  Systems 
Group  for  the  Defense  Atomic  Support  Agency  under  Contract  DASA0.1-70-C-0135 , 
Blast  Wave  Bouixdary  Layers  and  Particle  Entrainment.  This  final  report 
consists  of  two  separate  parts.  Part  I  contains  the  results  of  the  boundary 
layer  study,  and  Part  II  describes  the  results  of  the  particle  entrainment 
investigation.  The  first  part  is  independent  of  the  second,  while  the 
second  part  occasionally  refers  to  the  first. 

V.  Quan  was  the  project  engineer,  and  R.M.  Traci  and  J.F.  Farr,  Jr. 
contributed  significantly  to  this  study.  J.E.  Melde  and  V.R.  Hyman  also 
contributed  to  the  numerical  results.  The  interest  and  support  of 
H.J.  Carpenter  and  C.W.  Busch  in  this  work  are  acknowledged.  The  DASA 
technical  monitor  was  C.B.  McFarland. 
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ABSTRACT 

The  laminar  boundary  layer  flowfields  generated  by  plane,  cylindrical, 
and  spherical  strong  blast  waves  are  studied.  Solutions  are  obtained 
using  two  different  methods.  In  the  first  method,  the  governing  equations 
are  reduced  to  ordinary  differential  equations  by  using  a  similarity  trans- 
1 annation  and  a  parametric  integral  technique.  In  the  second  method,  the 
.  ssumptions  of  quasi-steady  state  and  local  similarity  are  invoked.  Results 
obtained  using  each  method  are  presented.  A  method  of  solving  the  turbulent 
boundary  layer  equations  is  also  outlined  in  detail  in  this  report. 
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NOMKNGLATURF. 


A  =  constant  defined  by  Equation  (3-3) 

Cp  =  specific  heat  at  constant  pressure 
C  =  constant  related  to  explosion  strength 
F  =  function  given  by  Equation  (A-5) 

h  =  specific  enthalpy 

k  =  conductivity 

£  =  distance  in  n  at  which  the  velocity  reaches  (l-e'^)  of  the  local 

free-stream  value. 

=  distance  in  n  at  which  the  temperature  reaches  (l-e”^)  of  the  local 
free-stream  value. 

m  =  geometry  factor,  =  2/(3+o) 

p  =  pressure 

Pr  =  Prandtl  number,  =  pc  /k 

P 

q  =  heat  flux  to  the  wall 
w 

Q  =  time  scaling  factor  defined  by  Equation  (5-4) 

R  =  function  given  by  Equation  (A-4) 
t  =  time 

u  =  velocity  in  x-direction 

u^  =  shock  velocity 

V  =  velocity  in  y-direction 

V  =  function  defined  by  Equation  (A-6) 

W  =  time  scaling  factor  defined  by  Equation  (5-7) 

X  =  distance  along  surface  from  explosion  point 

x^  =  distance  between  shock  and  explosion  point 
y  =  distance  normal  to  surface 
y^  =  boundary  layer  velocity  thickness 
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NOMENCLATURE  (Continued) 

=  boundary  layer  temperature  thickness 
a  =  geometry  factor,  =  (m-l)/m 

y  =  specific  heat  ratio 

6  =  distance  in  n  at  which  local  free-stream  velocity  is  reached 

6^  =  distance  in  n  at  which  local  free-stream  temperature  is  reached 

n  =  transformed  coordinate  defined  by  Equation  (3-2b) 

M  =  viscosity 

C  =  transformed  coordinate  defined  by  Equation  (3-2a) 

p  =  density 

T  =  transformed  coordinate  defined  by  Equation  (3-2c) 

T  =  shear  stress  at  the  wall 
w 

0  =0  for  two-dimensional  boundary  layer  and  1  for  axisymmetric  boundary 

].ayer 

0  =  0,  1,  and  2  for  plane,  cylindrical,  and  spherical  shock,  respectively. 

(p  =  function  given  by  Equation  (A-3) 

ip  =  stream  function  defined  by  Equation  (3-1) 

Subscripts 

e  ==  free-stream  (at  edge  of  boundary  layer) 
w  ==  wall  condition 

n  =  derivative  with  respect  to  q 

5  =  derivative  with  respect  to  5 

“  =  ambient  condition 
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1.  INTRODUCTION 


The  occurence  of  an  intense  explosion  at  or  near  the  earth's  surface 
results  in  a  complex  fluid-surface  interaction.  The  main  characteristics 
of  such  an  explosion  are  indicated  in  Figure  1-1.  The  intensely  hot,  low 
density  region  of  air  and  vaporized  earth  form  what  is  known  as  the  fire¬ 
ball  shortly  after  burst.  Within  seconds  the  fireball  begins  to  rise 
due  to  its  buoyancy,  carrying  with  it  large  amounts  of  vaporized  and 
pulverized  surface  material.  As  it  rises  it  is  cooled  by  entraining 
ambient  air  and  by  radiation  to  the  surrounding  air  and  surface.  In  the 
meantime  the  explosion  is  communicated  to  the  surrounding  air  by  the 

advancing  shock  front.  A  more  extensive  description  of  the  explosion 

1* 

phenomena  is  given  by  Erode.  The  fluid-surface  interaction  can  result 

in  large  amounts  of  ground  material  being  lifted  into  the  air.  Several 

such  lofting  mechanisms,  including  vaporization,  crater  splashing,  thermal 

expansion,  elastic  rebound,  jetting  from  cracks,  and  aerodynamic  entrain- 

2 

ment,  have  been  identified  by  Trulio  and  others. 

Of  prime  interest  to  this  study  is  the  boundary  layer  produced  by  the 
spherically  expanding  shock  wave  flowfield.  In  order  to  assess  the  aero¬ 
dynamic  entrainment  of  surface  materials,  it  is  necessary  to  determine 
the  boundary  layer  shear  stress  on  the  surface.  Also,  boundary  layer 
properties  are  required  for  the  design  of  surface  and  subsurface  installa¬ 
tions  which  are  hardened  to  the  thermal  and  dynamic  environment  of  the 
blast  wave.  Antennas,  silo  closures,  and  air  entrainment  systems  are 
examples  of  such  facilities.  Knowing  the  boundary  layer  properties  such 
as  the  velocity  and  temperature  profiles,  one  can  calculate  the  aerodynamic 
forces  on  the  facilities.  Test  techniques  such  as  (proposed)  NEST 


A 
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(Nuclear  Explosion  Simulation  Technique)  require  information  about  the 
boundary  layer  characteristics  to  properly  evaluate  the  test  technique 
and  to  interpret  test  data.  To  adequately  determine  the  degradation  of 
underground  facilities  caused  by  the  blast  environment,  the  flow  properties 
of  air  which  enters  weapon  system  ducts  are  required.  In  addition,  the 
rate  of  heat  transfer  to  exposed  surfaces  can  be  estimated.  The  total 
heat  transfer  consists  of  the  heat  convected  from  the  hot  boundary  layer 
gases  and  the  heat  radiated  from  the  fireball.  Of  particular  importance 
to  the  calculation  of  the  heat  transfer  is  the  generation  of  airborne 
dust  by  entrainment  into  the  boundary  layer.  The  dust  cloud  so  generated 
may  block  a  significant  portion  of  the  fireball's  thermal  radiation  by 
absorption  or  dispersion,  thus  decreasing  the  total  heat  rate  to  exposed 
surfaces.  . 

In  the  present  study,  the  boundary  layer  behind  a  strong  shock  moving 
into  a  stationary  fluid  is  investigated.  The  solution  is  presented  as  a 
model  for  the  blast  wave  boundary  layer  caused  by  an  explosion  at  or  near 
the  earth’s  surface.  Although  the  primary  interest  here  is  on  spherical 
blast  waves,  the  solutions  are  formulated  such  that  they  are  applicable 
to  plane  and  cylindrical  as  well  as  spherical  shocks.  Since  the  boundary 
layer  flow  is  laminar  within  a  short  distance  from  the  shock  front  in 
which  the  shear  stress  and  heat  transfer  are  high  and  since  it  is  of 
interest  to  assess  the  character  of  the  boundary  layer  if  it  were  to 
remain  laminar  throughout  the  blast  wave,  the  major  effort  of  the  present 
study  is  devoted  to  laminar  flow.  However,  a  method  of  solving  the 
turbulent  boundary  layer  equations  has  been  conceived  and  formulated  and 
is  included  in  the  present  report. 

The  equations  governing  the  two-dimensional  transient  flow  in  a 
laminar  boundary  layer  are  presented  in  Section  2  of  this  report  and  are 
solved  using  two  distinct  methods.  The  first  method  involves  a  similarity 
transformation  coupled  with  a  parametric  integral  technique  and  is  describ¬ 
ed  in  detail  in  Section  3.  This  method  is  believed  to  yield  accurate 
boundary  layer  solutions  for  strong-shock  conditions  (i.e.,  when  the 
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ambient  pressure  is  negligible  compared  to  the  shock  overpressure).  This 
condition  is  satisfied  for  early  time  blast  waves  which  are  of  most 
practical  and  theoretical  interest.  The  second  solution  method  involves 
the  assumption  of  quasi-steady  state  and  local  similarity  and  is  discussed 
in  Section  4.  This  method  is  being  investigated  because  the  solution  can 
be  extended  to  non-ideal  blast  waves  including  weak  shocks.  Results 
obtained  using  the  present  analyses  are  illustrated  in  Section  5,  and  a 
discussion  of  the  present  work  and  conclusions  are  given  in  Section  6. 

A  complete  definition  of  the  flowfield  involves  determining  the 

inviscid  flow  and  the  boundary  layer  flow.  The  inviscid  flowfield, 

calculated  by  neglecting  the  presence  of  the  boundary  layer,  has  been 

analyzed  by  a  number  of  investigators  and  is  rather  well  defined.  For 

3 

ideal  gas  and  strong  shock  conditions,  Taylor  obtained  an  exact  numerical 

4 

solution  and  an  approximate  closed-form  solution  while  Sedov  obtained  an 
exact  closed-form  solution.  Sedov's  inviscid  solution  is  used  in  the 
present  study  and  is  presented  in  Appendix  A. 

The  laminar  boundary  layer  associated  with  strong  shocks  has  been 
investigated  by  Mirels  and  Hamman;^  however,  their  solution  is  limited 
to  the  region  very  close  to  the  shock.  The  method  of  analysis,  results, 
and  region  of  applicability  of  this  solution  is  described  in  Appendix  B. 

This  solution  provides  the  only  rigorous  source  of  comparison  to  the 
results  of  the  present  study,  and  such  a  comparison  is  discussed  in 
Sections  3  and  4.  Mirels^  and  Murdock^  have  applied  integral  techniques 
to  obtain  solutions  to  the  boundary  layer  behind  shocks  of  constant 
velocity.  Such  a  shock  has  a  uniform  flowfield  behind  it,  whereas  the 
flowfield  for  a  blast  wave  is  far  from  uniform.  Thus,  those  solutions 
do  not  apply  to  the  problem  of  interest  here.  To  adequately  define  the 
blast  wave  boundary  layer,  a  more  comprehensive  solution  than  those 
presently  available  in  the  literature  is  needed.  The  need  is  for  a 
boundary  layer  solution  valid  throughout  the  region  of  the  influence  of 
the  blast  wave,  from  the  explosive  source  to  the  shock  front.  The  purpose 
of  this  study  is  to  present  such  a  solution  for  the  laminar  case.  Although 
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the  actual  blast  wave  boundary  layer  probably  becomes  turbulent  very  near 
the  shock  front  and  boundary  layer  theory  probably  is  invalid  near  the 
explosion  point  where  the  density  is  extremely  low,  it  is  believed  that 
the  laminar  solution  is  a  necessary  first  step  in  the  description  of 
the  physical  problem.  For  turbulent  boundary  layers,  a  method  of  solution 

g 

as  derived  by  Quan  is  shown  in  Appendix  C.  The  actual  numerical  results 
for  turbulent  flow  will  be  computed  at  a  later  date. 
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FIGURE  1-1  BLAST  CONFIGURATION 
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2.  GOVERNING  EQUATIONS 

A  blast  wave  is  by  nature  unsteady  and  produces  a  high  temperature, 
high  velocity  flowfield.  Thus,  transient  and  compressibility  effects  must 
be  necessarily  be  considered  in  the  solution  of  its  associated  boundary 
layer.  The  unsteady,  compressible  laminar  boundary  layer  equations  form 
the  starting  point  for  the  solution  perfomed  in  this  report,  and  are  pre¬ 
sented  in  this  section.  A  schematic  of  the  boundary  layer  which  defines 
the  parameters  used  is  given  in  Figure  2-1.  Both  two-dimensional  (a=0)  and 
axisymmetric  (a=l)  boundary  layers  are  considered.  The  y  coordinate  is 
measured  perpendicular  to  the  surface  and  the  x  coordinate  is  measured 
along  the  surface  from  the  point  of  origin  of  the  blast.  A  point  explosion 
generates  an  axisymmetric  boundary  layer,  a  plane  explosion  generates  a 
two-dimensional  boundary  layer,  and  a  line  explosion  can  generate  either  an 
axisymmetric  or  a  two-dimensional  boundary  layer. 

5  9 

The  boundary  layer  equations  applicable  to  both  cases  are:  * 

Continuity 


1R  +  Jl.  ^(pux  )  9(pv)  _  f. 

3t  o  9x  9y 


(2-1) 


Momentum 


9  /  9u\ 

9y  Y  9yj 


(2-2) 


-  =  -L  1-1 +  /IhV 

Dt  Pr  dy\^dyj 


(2-3) 


where  u  and  v  are  flow  velocities  in  the  x  and  y  directions  respectively, 
t  is  time,  p  is  pressure,  p  is  density,  h  is  enthalpy,  y  is  viscosity, 
and  Pr  is  Prandtl  number.  The  subscript  e  denotes  free-stream  properties, 
the  subscript  “  denotes  ambient  fluid  properties,  and 
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_D_ 

Dt 


3x 


+ 


(2-4) 


Air  is  assumed  to  be  an  ideal  gas  with  an  equation  of  state  given  by 

p  =  [(Y-1)/Y]ph  (2-5) 

where  y  is  the  ratio  of  specific  heats.  In  Equation  (2-3)  it  is  assumed 
that  the  constant  pressure  specific  heat  (c^)  and  Prandtl  number  (yc^/k) 
are  constant  where  k  is  the  fluid  conductivity.  The  general  boundary 
conditions  for  the  governing  equations  are: 


u(x,0,t)  =v(x,0,t)  »  0,  h(x,0,t)  “  h  (x,t)  at  y  =  0 

w 

(2-6) 

u(x,“,t)  =  u^(x,t),  h(x,«,t)  =  h^(x,t)  at  y  =  » 

Where  the  subscript  w  refers  to  the  flow  properties  at  the  wall  (y=0) . 
The  location  and  velocity  of  a  strong  shock  are  given  by 

X  =  Ct™,  u  =  Cmt°  ^  (2-7) 

s  ’  s 


where  the  constant  C  is  related  to  the  energy  liberated  during  the 
explosion  (see  Appendix  A)  and  the  constant  m  is  given  by  m=2/ (3+o) 
where  a=0,  1,  and  2  for  a  plane,  cylindrical,  and  spherical  shock, 
respectively  (see  Figure  2-2) . 

If  a  similarity  variable  is  defined  by  ^=l-x/x  ,  the  inviscid  flow 

s 

properties  can  be  express  as 


P  /Pco  ""  P  /Poo  “ 

e  ”  s  e  00 


u 

e 


%  h^ 


Y 

Y-1 


u 


2  FOI 
I  R(0 


(2-8) 


where  the  functions  F,  R,  and  (j)  depend  on  o  and  y  as  well  as  on  5  and 
they  are  given  in  Appendix  A. 
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CYLINDRICAL  SHOCK  SPHERICAL  SHOCK 

(a=a=1)  (7=2,  a=l) 


FIGURE  2-2  PLANE,  CYLINDRICAL,  AND  SPHERICAL  MOVING 
SHOCKS  IN  VICINITY  OF  PLANE  WALLS 
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3.  INTEGRAL  SOLUTION 


The  equations  given  in  Section  2  which  describe  the  laminar  boundary 
layer  behind  a  shock  have  been  reduced  by  Mirels  and  Hamman^  to  two 
partial  differential  equations  involving  two  dependent  and  two  independent 
variables.  Equation  (2-1)  is  satisfied  by  a  function  such  that 


^  OD 

u  * - V 

a  3y 
px  •' 


px 


/ail'  ?  a  p 


(3-1) 


The  following  independent  variables  are  introduced: 


5  “  1  -  x/Ct 

.  ry 


m 


_e_ 


dy 


1/2 


T  =  t 


(3-2a) 

(3-2b) 

(3-2c) 


where  A  is  a  constant  defined  by 


A  «  2mF  /p 

Q  CQ  *  □ 


(3-3) 


The  independent  variables  are  assumed  to  have  the  following  form 

U/2, 


u  i 

e 

m 


•f(5,n) 


(3-4a) 


h/h^  =  pjp  -  g(5.n) 

and  the  boundary  conditions  are  taken  as 


(3-4b) 


f(£,0)  =  “  0,  f^(4,”)  =  1 


g(5,0)  =  0,  g(5,“)  -  1 

CONFIDENTIAL 


(3-5a) 

(3-5b) 
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Similar  to  the  treatment  in  Reference  5,  the  wall  temperature  has  been 
assumed  to  be  negligible  compared  to  the  free-stream  temperature. 


Using  the  equations  given  above,  Equations  (2-2)  and  (2-3)  reduce  to 


(1-4)^°  f  +  (n-<j)f)f  =  2C<  r(|)r  +  <#>fr 

nnn  ^  '  nn  |  C  K 


,  -  y(2a+a)n  : 


+  (l-C-^f^)l-l 


f  /J  1  >**  ) 


(3-6a) 


25  *  |f<{>^+<j)f^-  2  (2a+oi)ri 


1  e  +  • 


('!i  +  !i  -  M 

\  g  rF  R  / 


g  (3-6b) 


where  a=(ra-l)/m.  In  deriving  Equations  (3-6a)  and  (3-6b)  it  is  assumed 
that  viscosity  is  proportional  to  temperature.  Equations  (3-6a)  and 
(3-6b)  must  be  solved  simultaneously  for  f  and  g  in  terms  of  5  and  n. 
It  is  noted  that  the  independent  variable  i  has  been  eliminated. 


The  shear  and  heat  transfer  at  the  wall  can  be  calculated  from 


■  ■■(s). 


2F  p  u  (x  -x) 

0*^00  S  S 


P  u  f  (5,0) 
tx^y  e  e  nn 


(3-7a) 


Asy/w 

[n  1/2 

_ _ ] 

2F  p  u  (x  -x)  \  / 

0*^00  s  9  J  \  s  / 


Pe  "e 


(3.-7b) 
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The  wall  temperature,  which  has  been  neglected  in  computing  the  coeffi¬ 
cients  of  drag  and  heat  transfer,  can  be  determined  by  the  procedure  out¬ 
lined  by  Mlrels  and  Hamman.^ 

Mirels  and  Hamman  obtained  a  solution  by  expanding  the  variables  f 
and  g  in  power  series  of  5  (see  Appendix  B) .  Their  solution,  however, 
is  valid  only  within  a  short  distance  from  the  shock  front.  The  purpose 
of  employing  the  integral  technique  is  to  obtain  a  solution  that  may  be 
valid,  under  the  assumption  of  laminar  flow,  for  all  distances  from  the 
shock  front. 

3.1  POLYNOMIAL  PROFILES 

The  integral  method  of  obtaining  approximate  solutions  to  boundary 
layer  problems  has  been  explained  in  textbooks  (e.g.,  Reference  10).  In 
the  present  analysis,  f  and  g  are  represented  by  polynomials  in  n.  The 
coefficients  of  the  polynomials  are  chosen  to  satisfy  the  boundary  condi¬ 
tions.  Equations  (3-6a)  and  (3-6b)  are  integrated  from  ri=0  to  ri="«  This 
results  in  two  first-order  ordinary  differential  equations  describing 

6  and  6  as  functions  of  5  where  6  and  6  are  the  distances  of  n  at  which 
t  t 

f  -»■  1  and  g  1,  respectively,  ft  should  be  pointed  out  that  6  and  6^ 
are  not  the  physical  boundary  layer  thicknesses.  The  actual  boundary 
layer  thicknesse.s  are  the  distances  in  y  at  which  n  6  and  n  and 

can  be  found  from  Equation  (3-2b)  once  6  and  6^  are  determined. 

In  the  present  analysis,  fourth-degree  polynomials  are  chosen  to 
represent  f^  and  g 

=  f  (C.n)  =  a,  z  +  a,-,z^  +  a.z^  +  a.z^*  (3-8a) 

u  n  l  2.  3  4 

e 

=  g(C,n)  “  b^w  +  b^w^  +  ^4”^  (3-8b) 

e 

where  z=n/6  and  w=ri/<Sj..  The  a^,  b^,  6,  and  6^  are  functions  of  4* 
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Equations  (3"8a)  satisfies  the  condition  of  f  (C,,0).  Upon  integration, 
Equation  (3-8ft)  yields  an  additional  constant  which,  to  satisfy  the  condi¬ 
tion  of  f(£.,0)*’0,  turns  out  to  be  zero.  The  four  coefficients  a,  are 

J 

determined  by  the  four  conditions  of  f^(C,6)*l,  f^^(5,6)->0,  (t;  ,6  )=0 , 

and  ,0)“0.  The  second  and  third  conditions  are  consistent  with  the 

definition  of  6,  i.e.,  all  derivatives  of  f^  vanish  at  ri"5,  and  the  fourth 
condition  is  obtained  by  satisfying  Equation  (3-6a)  at  ri“0. 

Similarly,  Equation  (3-8b)  satisfies  the  condition  of  g(C,t))=0.  The 

four  coefficients  b^  are  determined  by  the  four  conditions  of  g(C,<5^)=l, 

g_(5i'S  )"0,  g  (t,d  )“0,  and  g  (C,0)”  -Pr (Y-l)R(f ^f ^  ,0) /yi’  where  the 

n  c  nn  t  nn  rin 

fourth  condition  is  derived  by  evaluating  Equation  (3-6b)  at  ri«0. 

The  resulting  expressions  for  f  and  g  are 


w  2  1  4  1  5, 

6(z  -  Y  2  +  y  2  ) 


0-9a) 


g  =  2w  -  2w^  +  -  -j  b2(w  -  3w^  +  3w^  -  w^)  (3-9b) 


where 


2Pr(Y-l)E(}>^A^ 


(3-10) 


2  -  ri/6((,),  w  n/6^(4),  A(0  «  6^/6 


(3-11) 


Integrating  Equations  (3-6a)  and  (3-6b)  from  n^O  to  ri"=“  and  using 
Equation  (3-'Ja)  and  (3-9b)  ,  one  obtains  two  coupled  first-order  ordinary 
differential  equations  governing  6  and  6^.  These  equations,  after  simpli¬ 
fication  by  a  vast  amount  of  algebra,  are  given  by 


_ 

dC  “  J_ 
10 


(3-12a) 
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M 


d6.  G  +  6 


♦N  -  ^  (1-S) 


25 


dC 


(3-12b) 


where 


74 


315  20 


(2o+a) 


10 


b  Al-i 

60  2  /  R4. 


(3-13a) 


n  =  -  2  _F_  + 3  7  .V 

^  6  F  \10  "  315 


(3-13b) 


^  ,  r,,  ,  .  <16  ^  2a+ci  ^ 

G  -  +  4>  ^  1 - ^  <S 


t(^  +  h  >>2) 

^  ^  ^2  -  i ") 


(3-13c) 


K  =  b 


(  2  M  i 

’\  6  dC  R  4’  F  y 


(3-13d) 


(l_5)2a 

^  Pr  F  6 

o  t 


-2  +  |i  -  -l^lb. 


3  356/2 


io  ^^2)" 


(fi6H 

(3-13e) 


and 


,,  2, 2  3  .4j  1  .5,1,  AM  £  A-1 

H  -  15  A  -  j_^o  ^  180  ■^3'  ^2^^  for  A  __  X 


(3-14a) 


"  =  “  ^  ^  +  ifo 


^  4  b,AN  for  A  L  1 


3  '  4  3  2 

140A  180A 


I  = 


2  a2  9  ,4  1  .5 

15  ^  140  ^  "  45  ^ 


i  '>2''(^ 


30  ^  140  ^  126 


I  =  - 


3  +  ^ 


10  15A  55^3  5^^ 


3  9 

+  — 


14A^  70A^  36A 


for  A  £  1 
(3-14-b) 
1 
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W  ^  A  J.  3  ,3  1.4 

^  15  ^  35  ^  ~  36  ^ 


for  A  <  1 


(3-14c) 


M  =  - 


J_  +  _2_ 

15a2 


^  +  1 


140A^  45A^ 


for  A  >  1 


M  _  1  A  1  a3  ^  1  .4 

^  30  ^  140  ^  504  ^ 


for  A  <  1 


(3-14d) 


^  20  _  .  2  ,  ,  .  3 


5 


15A"  lAA”  280A^  180A^ 


for  A  >  1 


Equations  (3-12a)  and  (3-12b)  along  with  Equations  (3-13a)  to  (3-14d) 
are  programmed  for  numerical  solution  using  a  standard  Runge-Kutta-Adams- 
Moulton  method.  The  initial  values  of  6  and  6^  are  obtained  by  requiring 
d6/d?  and  d6^/dC  to  be  finite  at  5=0.  The  free-stream  properties  (F,  R, 
and  (j))  used  are  those  given  in  Appendix  A;  F  ,  R  ,  and  (p  are  derivatives 
of  F,  R,  and  <fi,  respectively,  with  respect  to  5«  The  results  for 
Pr=0.72  and  y=1.4  have  been  obtained  for  all  four  geometries  illustrated 
in  Figure  2-2.  The  results  of  6  and  6^,  along  with  other  results  to  be 
obtained  in  the  following  section,  are  shown  in  Figures  3-1  through  3-4. 

It  should  be  noted  that  once  6  and  6^  are  determined,  any  physical  property 
in  the  laminar  boundary  layer  can  be  computed  using  the  appropriate 
equations . 

3.2  EXPONENTIAL  PROFILES 

In  Section  3.1,  the  velocity  and  temperature  profiles  are  represented 
by  fourth-degree  polynomials.  It  is  highly  useful  to  obtain  an  alternate 
solution,  which  is  probably  less  accurate  but  is  simpler  to  derive,  by 
employing  profiles  of  another  form.  A  simple  representation  that  satisfies 
the  boundary  conditions  is  an  exponential  profile.  If  the  results  to  be 
obtained  using  the  exponential  representation  resemble  those  obtained 
using  the  polynomial  representation,  then  one  may  infer  that  the  integral 
solution  is  probably  not  too  sensitive  to  the  particular  profile  representa¬ 
tion  employed  as  long  as  they  are  physically  realistic. 
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The  exponential  forms  of  the  velocity  and  temperature  profiles  are 


f  =  1  -  exp(-ri/£) 
h 


(3-15a) 


g  =  1  -  exp(-n/^j.)  (3-15 b) 

where  H  and  SL  are  functions  of  g  and  are  the  distances  in  n  at  which  the 
velocity  and  temperature,  respectively,  reach  (1-e  )  of  the  local  free- 

stream  values.  It  should  be  noted  that  !L  and  are  measures  of,  but  have 
values  less  than,  the  transformed  velocity  and  temperature  boundary  layer 
thicknesses,  respectively. 

Integrating  Equations  (3-6a)  and  (3-6b)  from  y=0  to  y=«>  and  using 
Equations  (3-15a)  and  (3-15b) ,  one  obtains  two  coupled  first-order 
differential  equations  governing  Si  and  These  equations  are  given  by 


(1 


5  - 


2a4-a 

2 


where  Q=Si^/Si. 

Equations  (3-16a)  and  (3-16b)  are  programmed  for  numerical  solution 
using  the  Runge-Kutta  method.  The  initial  values  of  i  and  are  obtained 
by  requiring  dz/dg  and  dZ^/dg  to  be  finite  at  5=0.  The  free-stream 
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properties  used  are  those  given  in  Appendix  A.  The  results  for  Pr=0.72 
and  y=1«4  have  been  obtained  for  all  four  geometries  illustrated  in 
Figure  2-2.  The  results  of  and  which  are  sufficient  to  detenaine 
all  boundary  layer  properties,  are  shown  in  Figures  3-1  through  3-4. 

3.3  DISCUSSION  OF  INTEGRAL  SOLUTION 

The  integral  solution  for  laminar  boundary  layer  has  been  obtained 
using  two  distinct  representations  for  the  flowfield.  The  exponential 
profile  representation  is  simpler  for  computation  than  the  polynomial 
representation;  but  based  on  experience  with  simple  problems  such  as  flat 
plate  flow,  the  former  is  also  expected  to  be  less  accurate.  The  velocity 
and  temperature  profiles  obtained  using  the  two  methods  are  compared  in 
Figures  3-5  to  3-8,  respectively,  for  5=0,  0.02,  0.1,  and  0.5  and  spherical 
shock.  Mirels  and  Hamman's  solution,^  which  is  exact  as  5-^0,  is  also 
included  for  comparison  in  Figures  3-5  and  3-6.  Beyond  5*0.02,  however, 
the  accuracy  of  their  solution  is  questionable  (see  Appendix  B)  and  is 
therefore  not  shown  for  comparison  in  Figures  3-7  and  3-8.  The  factors 
that  determine  the  shear  stress  and  heat  transfer  at  the  wall,  f^^(5,0) 
and  gj^(5,0),  for  spherical  shock  are  plotted  in  Figure  3-9.  They  are 
magnified  and  compared  to  Mirels  and  Hamman's  solution  in  Figure  3-10 
for  small  values  of  5.  Also,  f^^(0,0),  g^(0,0),  [df^^ (5 >0) /d5]^  and 
[dg^(5 ,0) /d5] ^  for  plane,  cylindrical,  and  spherical  shocks  are  tabulated 
in  Table  3-1  below;  the  results  of  Mirels  and  Hamman's  are  Included  for 
comparison.  The  values  of  f^^(0,0)  and  g^(0)  are  independent  of  o  and  a. 

From  Figures  3-5  through  3-10  and  Table  3-1,  it  is  seen  that  the 
exponential  profile  solution  agrees  closely  with  the  polynomial  profile 
solution.  It  is  also  seen  that  the  polynomial  solution  yields  better 
overall  agreement  with  Mirels  and  Hamman's  solution  in  the  region  where 
their  solution  is  applicable. 
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TABLE  3-1 

COMPARISON  OF  SHEAR  STRESS  AND  HEAT  TRANSFER  FACTORS 
NEAR  SHOCK  FRONT  FOR  Pr  =  0.72  AND  y  =  1.4 


Mirels 

and 

Hamman 


Integral- 

Polynomial 


Integral- 

Exponential 


f  (0,0) 
nn 

g^(0,0) 


0.66141 

0.89693 


0.63579 

0.88779 


0.76376 

0.89032 


?,0) 


,0) 


0 

0 

1 

1 


-1/2 

-1 

-1 

-3/2 


1.1819 

1.7201 

2.5051 

3.0432 


1.1167 

1.6139 

2.4036 

2.9008 


-1/2 

-1 

-1 

-3/2 


-1.2891 

-2.7270 

-1.6291 

-3.0671 


-1.3824 

-2.8319 

-1.7548 

-3.2044 


1.6487 

2.2971 

3.3336 

3.9820 


-1.4124 

-3.1753 

-2.0089 

-3.7718 
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FIGURE  3-5  VELOCITY  AND  TEMPERATURE  PROFILES  AT  SHOCK  FRONT 
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FIGURE  3-6  VELOCITY  PND  TEMPERATURE  PROFILES  AT  1=0.02 
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4.  QUASI- STEADY  SOLUTION 


The  second  method  of  solving  the  blast  wave  laminar  boundary  layer 
is  referred  to  as  the  quasi-steady  solution.  Two  solution  procedures 
performed  in  a  quasi-steady  sense,  along  with  rationale  for  attempting 
them  are  discussed  in  this  section.  The  two  procedures,  though  independent 
in  method  and  application,  are  complementary  in  that  they  represent  vary¬ 
ing  degrees  of  approximation  to  the  actual  physical  problem.  The  first 
quasi-steady  procedure,  referred  to  herein  as  the  "steady  and  locally 
square  wave"  procedure,  is  described  in  Section  4.1.  This  procedure  uses 
available  solutions  for  the  boundary  layer  behind  a  square  wave  shock 
and  applies  them  in  a  locally  similar  sense  to  the  decaying,  spherically 
expanding  blast  wave.  The  "quasi-steady  and  locally  similar"  shock 
boundary  layer  solution  is  discussed  in  Section  4.2.  This  solution 
represents  a  more  complete  description  of  the  local  free-stream  flow 
properties  and  therefore  provides  a  more  realistic  test  of  the  appropriate¬ 
ness  of  the  local  similarity  assumption.  Both  solution  schemes  are 
developed  for  the  Taylor-Sedov  self-similar  blast  wave  flow  and  are  compared 
near  the  shock  front  to  Mirels  and  Hamman's  boundary  layer  solution 
(Appendix  B) .  A  discussion  of  the  applicability  of  both  solutions  to  the 
self-similar  blast  wave  and  possible  extension  to  non-ideal  blast  waves 
is  presented  in  Section  4.3. 

4.1  STEADY  AND  LOCALLY  SQUARE-WAVE  PROCEDURE 

The  "steady  and  locally  square  wave"  solution  procedure  was  an  obvious 
initial  step  in  examining  quasi-steady,  locally  similar  type  solutions 
to  the  blast  wave  boundary  layer  problem.  It  was  originally  performed 
with  the  intention  of  determining  the  accuracy  of  a  very  crude  local 
similarity  assumption,  based  on  known  shock  boundary  layer  solutions. 

The  comparison  of  the  resulting  solution  to  the  other  solutions  presented 
in  this  document  was  such  that  It  was  deemed  proper  to  briefly  discuss  the 
solution  method  and  results.  The  solution  is  discussed  in  the  following 
paragraphs  with  the  intention  of  suggesting  a  simple  procedure  for  getting 
rough  estimates  for  the  general  blast  wave  boundary  layer. 
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The  complete  solution  for  a  square  wave  shock  boundary  layer  is  given 
by  Mirels  in  References  6  and  11.  The  solution  was  obtained  by  solving 
the  boundary  layer  equations  in  a  steady,  shock  fixed  coordinate  system. 
The  solution  was  effected  by  reducing  the  equations  to  non-linear  ordinary 
differential  equations  with  a  similarity  transformation  based  on  the 
similarity  parameter: 

q  T2(x^-x)„^ 


-2-  dy 

P_ 


(4-1) 


The  subscript  q  refers  to  the  square  wave  boundary  layer  solution. 

A  numerical  solution  of  the  resulting  two  point  boundary  value  problem 

was  then  obtained  for  various  strength  shocks  as  defined  by  the  ratio  of 

the  shock  velocity  to  fluid  velocity  in  the  shock  fixed  coordinates  or 

u  / (u  -u  ).  These  solutions  are  tabulated  in  References  6  and  11.  The 
s  s  e 

procedure  for  applying  these  square-wave  solutions  to  the  boundary  layer 

for  a  transient,  spacially  decaying  blast  wave  is  to  use  the  tabulated 

solution  for  the  local,  instantaneous  value  of  the  velocity  parameter 

u  /(u  -u  ).  This  parameter  is  calculated  from  axiy  blast  wave  solution 
S  S  6 

(Taylor-Sedov  or  numerical)  using  the  instantaneous  shock  velocity  at 
the  time  of  interest  and  using  the  free-stream  velocity  (u^)  at  the 
position  behind  the  shock  of  interest.  In  this  simplified  sense,  then, 
the  solution  is  both  quasi-steady  and  locally  similar;  that  is  the  boundary 
layer  solution  is  assumed  to  depend  only  upon  the  instantaneous  shock  front 
and  local  free-stream  properties  (u  ,  x  ,  u  ,  T  ,  p  ) . 

S  S  6  3  6 

The  relevant  boundary  layer  parameters  from  Mirels  solution  are 
written  below  in  terms  of  the  variables  used  in  the  present  report.  The 
boundary  layer  velocity  and  enthalpy  are  given  by: 


u  -u 
s 


u  -u 
s  e 


q  q 


(4-2) 


where  the  prime  denotes  differentiation  with  respect  to  the  independent 
variable  h  and: 

q 
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h-  1  +  izi 


u 

— - —  M  t 

u  -u  e 
s  e  J  L 


r(n  )  -  r(0)s(n  ) 

q  q 


']  <t  ■  *) 


=<v 


2(u  -u  ) 

—  2.  s  e 

2u  +  (y-1)u 
s  e 


(4-3) 


The  wall  shear  stress  (t  )  and  heat  transfer  (q  )  are  given  by: 

w  ^ 


V(u  -  u  ) 

Y  '^7 — f  '  (0) 
2v  (x  -x)  q 
w  s 


(4-4) 


q  =  k  V^— 7— ^  s’ (0)<T  + 

w f  2v  (x  -x)  i  e 

*  w  s  { 


1  r  u  ~i' 

_ e _ - 

2  (u  -u  )  e 
s  e 


r(0)T  -  T 

e  w 


(4-5) 


The  boundary  layer  thickness  is: 


^  ^s~^e  >  3 

y«T2v^(x^-x)  ■  T„  ("S  10 


X  n*  +  f  ’ 
2  6  q 


where 


(4-6) 


V315CU  -  u  )‘ 

_ s  e 

189u  -  74u 

s  e 


and  where  y„  is  the  physical  momentum  boundary  layer  thickness.  The 


boundary  layer  thickness  is  derived  in  Reference  11  by  using  an  approximate 
integral  technique  which  compares  with  the  exact  numerical  technique 


within  5%, 


Curves  for  the  wall  shear  stress  parameter  [f^(0)]  and  heat  transfer 
parameter  [s'(0)]  and  recovery  factor  tt(0)]  are  plotted  as  functions  of 


the  free-stream  velocity  parameter  in  Figure  4-1  below.  Data  presented 
in  Reference  11  was  used  in  these  plots. 


For  purposes  of  comparison  with  Mirels  and  Hamman’s  solution,  the 
solution  procedure  outlined  above  was  completed  for  the  Taylor-Sedov  blast 
wave  flow  given  in  Appendix  A.  The  assumptions  made  in  the  square  wave 
boundary  layer  solution  are  consistent  with  the  assumptions  of  Mirels  and 
Hamman  for  the  self  similar  blast  wave  boundary  layer  (i.e.,  constant  Pr 
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and  Cp  and  y  linearly  proportioned  to  T)  so  that  these  assumptions  will 
not  affect  the  validity  of  the  comparison.  To  compare  the  two  solutions 
requires  a  transformation  from  the  square  wave  similarity  variable 
[Equation  (4-1)]  and  the  Mirels  and  Hamman  similarity  parameter.  Equations 
relating  the  relevant  parameters  are  given  below.  The  similarity 
variables  are  related  by: 


n 


q 


(1-0° 


2 

(y+1) 


F(0 


1/2 


n 


(4-7) 


The  non-dimensional  velocity  [f’(ri)]  and  non-dimensional  temperatures 
[8(n)]  are  related  through  Equations  (4-8)  and  (4-9): 

u  u  -  u 

f  ’  (n)  =  -^  -  — - -  f  ’  (n  )  (4-8) 

u  u  q 

e  e  ^ 

g(n)  =  1  +  ^  ~ 

(4-9) 


It  is  reiterated  that  the  prime  denotes  differentiation  with  respect 
to  the  independent  variable  (n  or  n^) .  The  parameters  which  determine 
the  shear  and  heat  transfer  are  related  through  Equations  (4-10)  and 
(4-11). 

(u  -u  )  /8n  \ 

- (4-10) 

g'  (n)  .  -  1  -  ^ 

where  9ri  /9ri  is  derived  from  Equation  (4-7). 

q 

A  meaningful  comparison  of  the  two  solutions  can  be  made  by  comparing 
the  values  of  f''(0)  and  g'(0)  as  calculated  by  Mirels  and  Hamman  and  as 
calculated  from  the  steady  and  locally  square  wave  procedures.  This 
essentially  compares  the  wall  shear  stress  and  heat  transfer.  The  values 
of  f''(0)  and  g'(0)  are  plotted  as  functions  of  the  non-dimensional 
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distance  behind  the  shock  front  in  Figures  4-2  and  4-3.  The  results  of 
the  steady  and  locally  similar  procedure  are  compared  to  Mirels  and  Hanman's 
two  term  expansion  (with  respect  to  £;)  very  near  the  shock  front  in 
Figure  4-2.  As  shown  in  that  figure  the  solutions  are  identical  at  the 
shock  front.  This  is  not  very  surprising  since  a  square  wave  is  the 
"zeroth"  order  approximation  to  the  decaying  self  similar  blast  wave  flow- 
field  which  is  changing  very  rapidly  with  C  near  the  shock  front.  The 
reasons  for  this  will  become  more  apparent  when  the  terms  in  the  boundary 
layer  formulation  which  are  neglected  in  the  local  square  wave  procedure 
are  examined  in  Section  4.2.  As  will  be  shown  in  Section  4.3,  the  accuracy 
of  the  solution  improves  with  distance  behind  the  shock  as  compared  to 
the  integral  solution  and  the  quasi-steady  and  locally  similar  solution. 

Some  reasons  for  this  will  be  discussed  in  that  section. 

4.2  QUASI-STEADY  AND  LOCALLY  SIMILAR  PROCEDURE 

The  "quasi-steady,  locally  similar"  blast  boundary  layer  solution  is 
presented  in  this  section.  The  spirit  of  the  solution  is  to  extend  Mirels 
and  Hamman's  solution^  to  all  regions  of  the  blast  wave  and  to  allow  for 
the  use  of  non-ideal  blast  waves  (weak  shocks,  precursors,  heights  of  burst) 
in  the  boundary  layer  solution.  This  section  includes  a  description  of 
the  solution  method  and  presents  results  using  the  Taylor-Sedov  self¬ 
similar  invlscld  blast  flowfieid.  Comment  on  the  accuracy  of  the  assump¬ 
tions  is  reserved  for  Section  4.3. 

The  basic  boundary  layer  equations  are  given  in  Section  2.  The 
first  step  in  the  solution  is  to  transform  the  continuity,  momentum  and 
energy  equations  [Equations  (2-1), (2-2),  and  (2-3),  respectively]  to  a 
shock  fixed  coordinate  system  and  to  concurrently  apply  a  similarity 
transformation  for  the  y  coordinate.  The  transformation  equations  used 
are  identical  with  those  of  Mirels  and  Hamman  and  are  given  by  Equation 
(3-2)  of  this  report.  The  effect  of  the  transformation  is  to  reduce  the 
effect  of  the  time  varying  terms,  thereby  making  the  quasi-steady  assump¬ 
tion  valid,  and  to  stretch  the  y  coordinate  to  account  for  compressibility 
effects.  As  was  pointed  out  by  Mirels  and  Hamman,  time  is  eliminated 
completely  by  their  similarity  transformation  for  the  self-similar  blast 
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wave  flowfield.  Thus  for  this  "ideal"  blast  wave  the  quasi-steady  assump¬ 
tion  is  superfluous.  This,  however,  is  not  the  case  for  a  "non-ideal" 
blast  wave. 

The  flow  velocity  components  (u,  v)  are  assumed  to  depend  upon  some 
scalar  function  '|'(5,h,T)  in  such  a  way  that  continuity  is  identically 
satisfied  as  given  by  Equations  (3-1)  earlier.  The  non-dimensional 
boundary  layer  velocity  and  enthalpy  are  thus  given  by: 


^  -  f(§,n)  ,  ^  =  g(5,n) 


(4-12) 


The  momentum  and  energy  equations  [Equations  (2-2)  and  (2-3)]  are  now 
transformed  as  described  above  into  Equations  (4-13)  and  (4-14)  below: 


2a.,  /  \  2a>,  v 

X  (1-c)  (py)  (1-5)  (py)„ 

_® _  f  a.  _  —5 - 0.  f 

.  2  2m(a+l)-l  nnn  *  2  2m(a+l)-l  nq 

Ap«T  Ap^T 

J.  9p  .  ,  u  . 

- S - ^  g  -.j(m(a+l)-l/2]  ^  ^  (1-5))  nf 

p  u  X  35  r  T  2x  f  nn 

e  e  8  s 

+  —  tA  -  [u  f  -  u  (1-5)1  f  r 
x^  5  nn  x„  e  n  s  n5 

s  S'-  J 


X 


9u 

_ e 

35 


(4-13) 
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2a.  ,-,2a,  . 

^  Xg  (1-C)  (pu) 

Pr  2  2m(a+l)-l  ®nn 

Apt 


1 

“  A  p2,2in(a+l)-l 


Pr 


1  1  ^Pe\ 

w 

“  /l  r  \ 


u  c  u  e  ru 

+  —  g  -  —  —  f  -(1- 

X  S  ri  X  u  n 

S  S  L  s 


-  |[m(a+l)-l/2]  f  ^  (l-0}n 


u  i  /ah 

e  / _ e 

■  h  x^  \  a  C 

e  s  \ 


(4-14) 


where  the  boundary  conditions  are  taken  as: 


f(c,0)  =  f  (5,0)  =  0  ,  f  (5,“)  =  1 

n  n 

(4-15) 

h 

8(5,0)  =  ~  ,  g(C,®)  “  1 


The  constant  A  in  the  above  equations  is  defined  by  Mirels  and  Hamman  and 

given  by  Equation  (3-3) .  For  a  general  inviscid  flowfield  the  constant  C 

is  defined  as  x  / t°*  so  that  A  is  represented  in  a  slightly  different  form 
s 

as : 


2mF  p 

o  s _ ^ 

^2m(a+l) 


(4-16) 
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It  is  important  to  note  that  the  governing  differential  equations 
(4-13)  and  (4-14)  reduce  identically  to  Mirels  and  Hamman's  equations 
(3-6a)  and(3-6b)  for  the  self-similar  blast  wave.  For  that  ideal  case, 
the  coefficients  in  Equations  (4-13)  and  (4-14)  become  independent  of  time. 
For  a  general  blast  wave  flowfield  these  coefficients  are  functions  of 
time  so  that  the  quasi-steady  assumption  is  required.  By  this  assumption 
then  the  solution  (f,  g)  is  a  function  of  time  but  only  in  the  sense 
that  the  coefficients  in  Equations  (4-13)  and  (4-14)  vary  with  time.  The 
quasi-steady  assumption,  in  essence,  involves  neglecting  all  terms  which 
involve  time  derivatives  of  the  functions  f  and  g.  As  noted  above  this 
is  justified  apriorl  for  the  Taylor-Sedov  blast  wave  by  the  elimination  of 
time  dependence  from  the  differential  equations.  The  assumption  essentially 
freezes  the  boundary  layer  at  a  given  Instant  of  time  and  the  solution  is 
assumed  to  depend  only  upon  the  instantaneous  flowfield.  Physically  the 
assumption  is  accurate  if  the  rate  of  diffusion  of  momentum  and  heat  through 
the  transformed  boundary  layer  is  greater  than  the  time  rate  of  change  of 
the  external  flow  properties.  Some  consequences  of  this  assumption  are 
discussed  in  Section  4.3. 

The  problem,  as  it  now  stands,  requires  the  quite  foxrmidable  solution 
of  a  pair  of  simultaneous  non-linear  partial  differential  equations. 

To  reduce  these  to  a  tractable  form  the  assumption  of  local  similarity 
is  made.  This  assumption  is  put  forth  in  the  same  cavalier  spirit  as 
the  quasi-steady  assumption.  Whereas  in  the  quasi-steady  assumption  time 
derivatives  of  f  and  g  were  neglected,  local  similarity  dictates  that 
space  derivatives  with  respect  to  5  of  f  and  g  are  neglected.  This  assump¬ 
tion  is  advanced  with  somewhat  weaker  justification  than  the  quasi-steady 
assumption.  Mirels  and  Hamman  have  shown  that  for  the  ideal  blast  wave, 
the  proper  similarity  transformation  makes  the  quasi-steady  assumption 
redundant.  Local  similarity,  however,  requires  neglecting  the  underlined 
terms  in  Equations  (4-13)  and  (4-14).  The  feeling  is  that  away  from  the 
shock  front  the  terms  neglected  by  local  similarity  become  smaller  relative 
to  the  other  terms  in  Equations  (4-13)  and  (4-14) .  The  fair  agreement 
with  the  integral  solution  as  seen  in  Section  4.3  is  advanced  as  justifica¬ 
tion  for  the  local  similarity  assumption. 
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In  essence  the  assumption  of  local  similarity  requires  that  all 
boundary  layer  flow  parameters  (u/u^,  h/h^,  pp)  be  a  function  only  of  the 
similarity  variable  n  at  any  given  time  (t)  or  distance  behind  the  shock 
(C).  This  assumption  reduces  Equations  (4-13)  and  (4-14)  to  two  simul¬ 
taneous  non-linear  ordinary  differential  equations  with  split  boundary 
conditions.  The  resulting  two-point  boundary  value  problem  can  be  solved 
by  standard  numerical  methods.  The  method  used  is  known  as  a  shooting 
technique  wheteby  guesses  are  made  for  the  higher  order  derivatives  at 
the  wall  [f^^(C,0),  g^(C,0)].  Equations  (4-13)  and  (4-] 4)  are  then 
integrated  numerically  to  some  large  value  of  n.  Guesses  are  then  made 
for  f^^(4,0),  g^((;iO)  until  the  resulting  solution  gives  values  of 
f^(5.°“)  and  g(^,“)  which  are  equal  to  1  within  some  arbitrary  small  error. 
This  can  be  carried  out  by  an  iteration  scheme  which  is  found  to  converge 
quite  rapidly  and  efficiently  to  the  proper  solution  This  process  is 
repeated  for  the  values  of  C“l-x/x^  of  interest.  This  can  be  performed 
for  the  blast  wave  flowfleld  at  any  instant  of  time  thereby  determining 
the  blast  boundary  layer  as  a  function  of  position  (x,y)  and  time  t. 


Once  the  functions  t  and  g  (and  their  derivatives)  are  determined  for 
the  values  of  K  and  x  of  interest  the  problem  is  essentially  solved.  It 
remains  to  calculate  the  physical  parameters  of  interest  and  re-transform 
the  independent  variables  C.  n  back  into  physical  coordinates  (x,y). 

This  is  done  using  the  following  equations: 


Xgd-C) 


y  = 


At 


2m(a+l)- 


1/2 


rn 


g(?,n)dn  (4-17) 


The  velocity  and  enthalpy  become: 

u  =  u^(x,y,t)f^[5(x,t)  ,  n(x,y,t)] 
h  =  h^(x,y,t)g[C(x,t)  ,  n(x,y,t)] 


(4-18) 
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Parameters  of  particular  interest  are  the  wall  shear  stress  and  heat 
transfer.  These  are  given  by  Equations  (4-19)  and  (4-20),  respectively, 
r  -  1/2 

/  X  V 

'w'  - -  (^1 


V  Pr  . 


o  V  2  m-1,  , 


2 

(fj 


(c,0)  (4-20) 


Equations  (4-13)  and  (4-14)  were  solved  using  the  procedure  outlined 
above  for  the  Taylor-Sedov  blast  wave  flowfield  and  by  assuming  Pr,  c^ , 
and  m/T  to  be  constant.  The  results  for  the  shear  stress  and  heat  transfer 
parameters  are  given  near  the  shock  front  in  Figure  4-4  and  for  the  complete 
flowfield  in  Figure  4-5.  As  expected  these  parameters  compare  more 
favorably  with  Mirels  and  Hamman's  solution  near  the  shock  front  than  the 
steady  and  locally  square-wave  solution.  It  is  noted  that  the  solution  is 
identical  with  Mirels  at  the  shock  front  and  that  very  near  shock  front 
the  trend  in  both  f^^  and  g^  compare  favorably.  Boundary  layer  velocity 
and  enthalpy  profiles  at  the  shock  front  and  at  ^-0.5  are  given  in 
Figures  4-6  and  4-7.  Of  note  in  these  figures  is  the  indicated  decrease 
in  transformed  boundary  layer  thickness  with  distance  behind  the  shock 
front . 


4.3  DISCUSSION  OF  QUASI-STEADY  SOLUTIONS 

Figures  4-8  through  4-11  provide  a  comparison  of  the  wall  shear 
stress  and  heat  transfer  parameters  as  calculated  by  each  of  the  flow 
calculational  schemes  discussed  in  this  report.  Figures  4-8  and  4-9 
compare  the  solutions  very  near  the  shock  front  in  which  region  Mirels 
and  Hamman's  original  boundary  layer  solution  is  applicable,  and  Figures 
4-10  and  4-11  compare  the  solutions  presented  in  this  report  over  the 
entire  blast  wave.  Figures  4-8  and  4-9  show  that  all  the  solutions 
compare  well  near  the  shock  front  with  the  integral  solution  showing  the 
best  agreement  with  Mirels  and  Hamman's  solution  with  respect  to  variation 
of  f  and  g  with  distance  from  the  shock  front  (C) .  This 

nn  n 

trend  with  5  is  indicative  of  the  inadequacy  of  the  assumptions  of  local 
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similarity  of  both  quasi-steady  solutions.  This  point  is  further  veri¬ 
fied  by  the  fact  that  the  quasi-steady  and  locally  similar  solution  shows 
a  slightly  better  trend  with  respect  to  ^  than  the  steady  square  wave 
solution.  This  is  attributed  to  the  better  approximation  of  variations 
with  4  in  the  quasi-steady  and  locally  similar  formulation.  The  rather 
good  comparison  of  both  quasi-steady  solutions  with  the  integral  solution 
away  from  the  shock  front  is  indication  that  the  assumption  of  local 
similarity  improves  with  distance  from  the  shock  front.  This  is  to  be 
expected  since  the  blast  wave  flowfield  (and  boundary  layer  flow)  varies 
less  rapidly  with  respect  to  4  as  4  increases  away  from  the  shock  front. 
It  is  thus  expected  that  the  terms  neglected  by  local  similarity  (f^ 
and  f^^)  will  decrease  in  significance  (relative  to  the  other  terms  in 
the  boundary  layer  equation)  as  4  increases. 

It  is  difficult  to  justify  or  to  speculate  on  the  accuracy  of  the 
steady  and  locally  square  wave  solution.  Suffice  it  to  point  out  that 
the  above  figures  show  that  the  solution  compares  surprisingly  well  with 
the  more  rigorous  integral  solution  for  the  Taylor-Sedov  blast  wave. 

The  good  agreement  is  fortuitous.  The  solution's  utility  rests  on  the 
fact  that  it  provides  a  relatively  simple  procedure  for  estimating  the 
important  boundary  layer  flow  parameter  (shear  stress,  heat  transfer  and 
boundary  layer  thickness).  It  can  also  be  used  to  estimate  such  para¬ 
meters  for  the  boundary  layer  of  a  general  non-ideal  blast  wave.  The 
relatively  good  comparison  shown  here  with  the  more  rigorous  boundary 
layer  solutions  should  add  some  credibility  to  the  numbers  calculated 
in  such  a  manner. 

The  results  of  the  quasi-steady  and  locally  similar  solution  for  the 
Taylor-Sedov  blast  wave  are  quite  satisfying.  The  primary  advantage 
of  the  solution,  and  the  reason  for  performing  it,  is  that  it  is  not 
dependent  upon  a  self-similar  inviscid  flowfield,  i.e.,  the  assumption  of 
a  strong  shock  is  not  required.  This  opens  the  possibility  of  extending 
the  solution  to  include  non-ideal  effects  in  the  free-stream.  In 
particular,  a  numerical  solution  valid  in  the  weak  shock  regime  or  one 
which  includes  radiation  effects  (which  cause  a  shock  precusor)  could  be 
used  for  the  boundary  layer  analysis.  The  accuracy  of  the  quasi-steady 
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and  local  similarity  assumptions  for  such  a  flowfie.ld  can  only  be  con¬ 
jectured  upon.  The  quasi-steady  assumption  could  be  quite  good.  As 
pointed  out  earlier,  the  similarity  transformation  used  eliminates  time 
from  the  boundary  layer  equation.s  [Equations  (4-13)  and  (4-14)J  for 
any  power  law  shock  (x  =Ct™)  for  which  the  trailing  blast  wave  flow 
properties  have  a  similar  spatial  distribution  for  different  times.  This 
condition  at  least  partially  holds  for  a  high  explosive  or  nuclear  blast 
wave  for  even  late  time.  That  is,  even  for  weak  shocks,  the  spatial  dis¬ 
tribution  of  blast  wave  flow  properties  varies  little  in  general  shape. 
Thus  it  is  conceivable  that  even  for  such  a  case  the  quasi-steady  and 
locally  similar  procedure  could  give  good  results.  As  far  as  the 
applicability  of  the  local  similarity  assumptien  is  concerned,  it  can 
only  be  pointed  out  that  the  assumption  provides  quite  good  results  in 
comparison  to  the  integral  solution  for  the  Taylor-Sedov  flowfield. 

It  is  expected  that  this  would  carry  over  to  the  more  general  blast  wave. 

The  quasi-steady  and  locally  similar  solution  also  presents  the 
possibility  of  being  extended  to  include  turbulence  by  using  a  phenomeno¬ 
logical  theory  of  turbulence,  or  the  solution  could  be  numerically 
coupled  to  a  sami-empirical  theory  for  wall  turbulence.  The  solutions 
presented  in  this  report  assume  viscosity  is  a  linear  function  of 
temperature.  This  restriction  is  not;  required  in  the  quasi-steady  and 
locally  similar  solution.  In  any  practical  calculation  viscosity  would 
be  made  a  function  of  temperature  using  an  empirical  equation  like  the 
Sutherland  viscosity  relation.  This  would  be  consistent  with  local 
similarity  and  pu  would  be  some  function  of  the  temperature  ratio  g.  In 
summary,  this  solution  method  is  attractive  becau'e  of  the  straight¬ 
forward  manner  in  which  it  may  be  extended  to  include  effects  which  are 
difficult  to  handle  theoretically  but  which  add  considerably  to  the  ac¬ 
curacy  of  the  physical  model. 
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STEADY,  LOCALLY  SQUARE  WAVE  SOLUTION  WALL  SHEAR  STRESS  PAPA.METER 


STEADY,  LOCALLY  SQUARE  WAVE  SOLUTION  WALL  HEAT  TRANSFER  PARAMETER  [s'(0)]  AND  RECOVERY  FACTOR  [r(0)] 


FIGURE  4-4  WALL  SHEAR  AND  HEAT  TRANSFER  PARAMETERS  FOR  STRONG 
SPHERICALLY  EXPANDING  SHOCK  NEAR  SHOCK  FRONT  USING 
QUASI-STEADY,  LOCALLY  SIMILAR  SOLUTION 
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FIGURE  4-6  BOUNDARY  LAYER  VELOCITY  PROFILES 
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FIGURE  4-8  WALL  SHEAR  STRESS  PARAMETER  FOR  STRONG,  SPHERICALLY  EXPANDING  SHOCK  NEAR  SHOCK  FRONT 
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FIGURE  4-9  WALL  HEAT  TRANSFER  PARAMETER  FOR  STRONG,  SPHERICALLY  EXPANDING  SHOCK  NEAR  SHOCK  FRONT 


(Reverse  of  Page  is  Blank) 
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5.  ILLUSTRATION  OF  RESULTS 


The  solutions  to  the  blast  wave  laminar  boundary  layer  equations 
have  been  given  in  Sections  3  and  4.  Of  the  several  solutions  presented, 
the  most  rigorous  one  is  the  integral  solution  with  the  polynomial 
representation  for  velocity  and  temperature  profiles.  This  solution  is 
employed  to  illustrate  some  physical  properties  of  interest  in  this 
section. 

The  velocity  thickness  y^  and  the  temperature  thickness  y^  can  be 
obtained  by  setting  n=<S  and  respectively,  in  Equation  (3-2b).  The 

result  is 


t  \  /  (1-0°  ^ 


(5-1) 


^6  =  ^6 


y(io  -  60  ‘’2) 


(5-2) 


where 


1  _  J _ 1_  . 

A  10  60  2 


J; _ ^  _1^ 

ZA*^  5A^ 


_  1  .  /_A_  _  -L  +  ^ _ L_\ 

^  ^  \2A^  A^  4A^  5A^/ 


for  A  <  1 


(5-3) 


for  A  >  1 


Thus,  both  y^  and  y^  can  be  given  as  products  of  a  function  of  t  and  a 
function  of  C*  By  letting 


(5-4) 
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the  results  for  and  /Q(t)  are  plotted  in  Figures  5-1  and  5-2 

for  spherical  shock.  Since  A  is  proportional  to  the  scaling 

is  that  and  y^  are  proportional  to  C  which  varies  with  the  explosion 

energy  E  as  Also,  at  given  y.  and  y»  are  proportional  to 

in'" 0  5  0  4  0  0^ 

t  .  For  C=3800  ft/sec  '  (corresponding  roughly  to  a  1.0  megaton  yield*) 
and  t=0.1  and  1.0  sec,  y^  and  y^  are  shown  as  functions  of  x  in  Figure 
5-3.  It  is  seen  that  the  boundary  layer  thicknesses  are  quite  thin  for 
a  long  distance  behind  the  shock,  after  which  the  thicknesses  increase 
rapidly.  The  velocity  thickness  for  the  blast  wave  at  0.1  second  is 
compared  to  the  velocity  thickness  for  a  square  wave  (using  peak  uniform 
free-stream  properties)  in  Figure  5-4.  It  is  seen  that  large  differences 
exist. 

The  shear  stress  t  and  heat  transfer  rate  q  can  be  written  as 
w  ^w 


T  =  W(t) 
w 


a-o 


F+(|) 


(5-5) 


q„  =  w(t) (Cmt 
w 


/2  E  J 


(5-6) 


where 


(5-7) 


Thus  can  be  scaled  with  W(t),  and  q^  with  W(t)Cmt  or  'W(t)u  (t). 

From  the  expression  for  W(t),  it  is  seen  that  t  and  q  decrease  rapidly 

as  time  is  increased.  The  functions  'r^/W(t)  and  q^/W(t)u^(t)  are  plotted 

in  Figure  5-5  for  the  spherical  shock.  It  is  seen  that  near  the  shock, 

both  shear  stress  and  heat  transfer  rate  decrease  with  distance  from  the 

shock.  Some  distance  away,  however,  the  heat  transfer  rate  increases 

rapidly.  The  reason  for  this  is  that  the  free-stream  temperature 

increases  more  rapidly  than  the  heat  transfer  coefficient  decreases  in 
- 1 _ /-t 'i  Cir\r\  _ 0.4  j  ^  ~  ^  , 


this  regime.  For  C=3800  ft/sec  '  and 
shown  as  functions  of  x  in  Figure  5-6. 


and  t=0.1  and  1.0  sec,  x  and  q  are 

w  ^w 


*This  is  calculated  based  on  1.0  megaton  free  air  burst  which  corresponds 
roughly  to  1/2  megaton  surface  burst. 


CONFIDENTIAL 


j'  TwH  rj'.  r.-  ,v  r>  vh  "j(  nk  r r j.  *  >  r  >  fvv  r j..  fu. a >7- >  x  ^  k  ^  w  rx  >  x  KKT'X  > x  ?w-Xvv>x  t x  ivk  xxTUf  xvw 


CONFIDENTIAL 


16118-6002-R7-00 


The  velocity  and  temperature  distributions  shown  in  Section  3  are 
given  in  terms  of  n.  These  distributions  can  be  given  in  terms  of  y 
through  Equation  (3- 2b)  which  relates  n  to  y.  For  the  polynomial 
solution,  the  relation  becomes 


y  =  Q(t) 


.1/2 


(1-0 


1/2  R 


y  =  Q(t) 


-1/2 


(1-0 


1/2  R 


'^tllO 


60 


+  n 


6 

t 


for  n 

for  n  >  6^ 
—  t 


(5-8) 


The  results  of  the  velocity  and  temperature  profiles  as  functions  of 
y/Q(t)  and  5  are  shown  in  Figures  5-7  and  5-8  for  the  spherical  shock. 


I 

i 


The  following  and  last  illustration  does  not  require  the  results  of 
the  boundary  layer  study  but  may  be  of  interest.  It  is  difficult  to 
assess  the  transition  point  from  laminar  flow  to  turbulent  flow.  However, 
one  may  wish  to  estimate  the  Reynolds  number  variation.  For  this  purpose, 
one  may  define  a  Reynolds  number  Re  as 


u  (x  -x)  p  t 

Oft 


Re  = 


<j)R  i 


M  m 


p  t 

‘  CO 

y  m 


z(0 


(5-9) 


The  value  of  Re/ (p^t/p^m)  is  plotted  in  Figure  5-9  for  the  spherical  shock. 
It  is  seen  that  Re  first  increases  and  then  decreases  as  i  is  increased. 
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FIGURE  5-2  SCALED  BOUNDARY  LAYER  VELOCITY  AND  TEMPERATURE  THICKNESSES 
NEAR  SHOCK  FRONT 
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FIGURE  5-8  SCALED  TEMPERATURE  PROFILES 
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6.  DISCUSSION  AND  CONCLUSIONS  j 

j 

j 

The  blast  wave  laminar  boundary  layer  equations  have  been  solved  by  | 

four  procedures.  The  most  rigorous  solution  is  the  one  obtained  using  i 

i 

the  integral  technique  with  the  polynomial  representation  for  velocity 
and  temperature.  The  integral  solution  with  the  exponential  representation 

i 

is  simpler  to  derive  and  it  agrees  quite  will  with  the  integral  solution  i 

with  the  polynomial  representation.  The  steady  and  locally  square  wave  j 

method  allows  for  rapid  estimation  of  the  flow  properties;  however,  although  | 
it  yields  fair  agreement  with  the  more  rigorous  solutions,  it  lacks  theo-  ] 

retical  justification.  The  quasi-steady  and  locally  similar  method  is  a  more 
reasonable  approach  than  the  steady  and  locally  square  wave  method  and  the  i 

results  are  generally  more  accurate,  but  it  requires  more  extensive  computation., 

i 

i 

In  the  region  near  the  shock  front,  only  the  integral  solution  yields  i 

5  i 

close  agreement  with  a  previously  established  solution.  Since  the  integral 

! 

solution  is  obtained  using  a  rigorous  method,  this  solution  is  recommended  | 

for  general  use.  However,  if  weak-shock  and  non-ideal  gas  effects  must 
be  considered,  then  the  quasi-steady  method  can  be  used  to  obtain  rough  | 

estimates  of  the  flow  properties.  ; 

I 

i 

The  integral  solution  shows  that  the  boundary  layer  thicknesses 
remain  quite  thin  for  a  long  distance  behind  the  shock,  then  they  increase 

very  rapidly  with  the  distance.  This  is  in  contrast  to  a  flat  plate  or 

square-wave  solution  for  which  the  boundary  layer  thickness  increases 
with  distance  behind  the  shock  to  0.5  power.  The  integral  solution  also 
shows  the  temperature  thickness  to  be  much  greater  than  the  velocity 
thickness  except  for  a  short  distance  from  the  shock.  Near  the  shock 
front,  both  shear  stress  and  heat  transfer  decrease  as  the  distance 
behind  the  shock  is  increased.  However,  at  a  distance  behind  the  shock 
of  about  one-third  of  the  length  between  the  shock  and  the  explosion 
point,  the  heat  transfer  increases  very  rapidly.  This  is  due  to  the  fact 
that  the  free-stream  temperature  increases  at  a  much  higher  rate  than 
the  heat  transfer  coefficient  decreases  as  the  distance  behind  the  shock 
is  increased. 
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The  results  presented  in  Section  3  for  the  integral  solution  allow 
one  to  compute  any  physical  properties  for  laminar  flow  such  as  shear 
stress,  heat  transfer,  velocity  and  temperature  profiles,  etc.,  by  using 
the  appropriate  equations.  The  results  can  also  be  used  for  strong 
shocks  of  various  explosion  energies  by  using  proper  scaling.  Some 
illustrations  of  the  results  are  shown  in  Section  5  where  the  method  of 
scaling  is  also  indicated. 

The  flow  is  expected  to  be  laminar  only  for  a  short  distance  from 
the  shock.  It  then  becomes  turbulent  for  a  predominant  portion  of  the 
flowfield.  At  less  than  half  the  distance  between  the  shock  and  the 
explosion  point,  the  density  becomes  extremely  low  and  the  continuum 
boundary  layer  equations  may  become  invalid.  However,  the  boundary 
layer  properties  near  the  explosion  point  are  not  of  primary  interest. 
Hence  the  next  flow  regime  that  should  be  investigated  is  that  of 
turbulent  flow.  For  turbulent  boundary  layers,  a  method  of  solution 

g 

has  been  formulated  by  Quan  and  is  shown  in  Appendix  C.  It  is  recom¬ 
mended  that  calculations  be  carried  out  to  obtain  numerical  results. 
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APPENDIX  A 


BLAST  WAVE  INVISCID  FLOWFIELD 


The  inviscld  flowfleld  generated  by  an  intense  explosion  in  a  uniform 


atmosphere  has  been  studied  by  a  number  of  investigators.  In  a  classic 
3 

paper,  Taylor  presented  a  self-similar  solution  for  a  point  explosion. 


He  presented  an  exact  solution  in  numerical  form  and  an  approximate 

4 

solution  in  closed-form.  Sedov  independently  arrived  at  an  exact  solution 


and  generalized  it  to  a  line  explosion  (cylindrical  shock)  and  a  plane 
explosion  (plane  shock).  Sedov's  solution  is  given  in  analytical  form 
and  is  the  one  used  in  this  report. 


A  strong  shock  moves  according  to  a  power  law  in  time  t  of  the  form 


u  =  Cmt 
s 


(A-1) 


where  x  is  the  shock  position  and  u  is  the  shock  velocity.  The  constant 
s  s 


m  is  determined  by  shock  geometry  and  is  given  by  jn=2/(3+ar)  where  a=0,  1, 
and  2  for  plane,  cylindrical,  and  spherical  shocks,  respectively.  The 


constant  C  is  determined  by  the  explosion  strength  and  is  given  by 


C*=(E/p  )"*  where  p  is  the  ambient  density  and  E  is  a  certain  constant 


which  has  the  dimensions  of  and  is  proportional  to  the  energy  E^  (kinetic 


energy  and  heat  energy)  liberated  during  the  explosion,  l.e.,  E=aE^.  The 


constant  a,  which  depends  on  y  ai'd  a,  has  been  shown  in  graphical  form  by 


Sedov.  For  air  and  spherical  shock,  a=*i.i75.  it  should  be  pointed  out 


that  E^  is  not  the  total  energy  liberated  in  an  explosion  since  a  large 


part  of  the  energy  is  expended  in  radiation. 


The  inviscid  flow  properties  behind  the  shock  are  given  by 


p  =  p  u  F(0  j 


e 


(A-2) 


=  u  <(>(5)  . 
e  s 


h  =  - 

e  Y-1  P, 
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where  p,  p,  a,  and  h  are  the  pressure,  density,  velocity,  and  enthalpy, 

respectively.  The  subscript  e  refers  to  conditions  behind  the  shock, 

and  "  refers  to  the  undisturbed  atmosphere;  F,  R,  and  ()>  are  the  solution 

parameters  given  as  functions  of  the  non-dimensional  blast  wave  coordinate 

5=l-x/x  .  From  Sedov's  results,  one  may  write 
s 


(?) 


v(i-0 


(A- 3) 


R 


Y+1  / (v+2)y  V 

“3 

1±1  ( 

1  _  )^±2  V  ^ 

Y-1  V  2 

VJ 

[y-1  V 

2  )j 

(v+2)  (y-H) 


(v+2)(y+1)-2[2+v(y-1)] 
-1  2v  I- 


2  \ 

(v+2) (y+1) 

Y+i; 

4 

2+v 


(v+2)  (Y'bl) 


(v+2)(Y-l)-2t2+v(Y-l)] 
where  V  is  an  implicit  function  of  C  given  by 
5  =  1- 


.  2+vitzii  v) 


(A-4) 


2+v(y-1)  V 
2 


‘4“^“l 


(A-5) 


1  2 

n 

-] 

“  2+v 

xti 

^(v+2)y  y  ,) 

Y-1 

V  2  / 

-ct. 


(v+2) (y+1) 


(v+2)(y+1)-2[2+v(y-1)] 


1  - 


2+v(y-1)  y 
2 


) 


-a. 


(A-6) 


and 
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V  =1+0 


a,  = 


_  (v+2)y 


2+v(y-1) 

1-Y 

2(y-1)+v 

V 

3  "  2(y-1)+v 


2v(2-y) 

y(v+2)^ 


= 


Oi„  = 


a ,  = 


= 


a^(v+2) 
2  -Y 

2 


Y 


_0 


(A- 7) 


The  functions  (f> ,  R,  and  F  for  a=2  (spherical  shock)  and  Y'^I-A  are  shown  in 
Figure  A-1  which  illustrates  the  severe  nonuniformities  near  the  shock 
front . 


Only  the  similar  solution  as  given  in  this  appendix  is  used  in  this 
report,  since  the  integral  boundary  layer  solution  requires  the  similarity 
transformation.  It  should  be  mentioned,  howev . >  ,  that  non-similar  solu¬ 
tions  to  the  inviscid  flowfield  are  available  and  they  can  be  used  in  the 
quasi-steady  method  of  boundary  layer  solution.  The  non-similar  inviscid 
solutions  involve  the  numerical  solution  of  the  continuity,  momentum, 
and  energy  equations  which  take  into  account  the  real  gas  effects,  ambient 
pressure,  and  radiation  heat  transport.  In  the  strong  shock  phase, 
these  numerical  solutions  agree  closely  with  the  similar  solution. 
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APPENDIX  B 

MIRELS  AND  HAMMAN'S  BOUNDARY  LAYER  SOLUTION 

As  noted  In  Section  1,  Mirels  and  Haniman  performed  a  solution  to 
the  blast  wave  boundary  layer  which  is  valid  near  the  shock  front. ^ 

For  reference  purposes,  the  method  and  results  of  that  solution  are 
outlined  in  this  appendix. 

The  basic  transformed  equations  and  boundary  conditions  used  by 
Mirels  and  Hamman  are  given  in  Section  3.  These  equations  are  partial 
differential  equations  with  two  independent  variables  i  and  n.  Mirels 
and  Hamman' s  solution  method  is  to  expand  the  dependent  variables  f  and 
g  as  power  series  in  Substituting  f  and  g  into  the  differential 
equations  and  equating  coefficients  of  like  powers  of  ^  reduces  the 
equations  to  a  series  of  non-linear  ordinary  differential  equations 
which  can  be  integrated  using  a  shooting  technique  as  discussed  in 
Section  4.2.  The  solution  can  thus  be  continued  to  any  desired  order 
of  approximation  although  each  successive  approximation  becomes  increas¬ 
ingly  more  complicated. 

The  functions  f  and  g  are  expanded  in  the  form 

f(^,ri)  =  f  (n)  +  (a5)f  (n)  +  (a0^f,(n)  +  ... 
o  i  z 

(B-1) 

g(c,n)  =  +  (aC)gj^(ri)  +  (a0^g2(n)  +  ... 

Sedov's  solution  to  the  freestream  flowfield  parameters  is  expanded  in  a 
similar  power  series  to  the  same  degree  of  accuracy  as  follows: 

F(£:)  -  F  +  (a4)F^  +  (uC)^F„  +  ... 

o  1  2 

^(5)  =  +  (ae)c|)^  +  (aC)^(!)2 

R(5)  =  R  +  (aC)R.  +  (a^)X  +  ... 

o  1  2 
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The  constant  a  is  dependent  upon  the  exponent  In  the  power  law  shock 
relation  (Appendix  A)  with  the  fom  a=(m-l)/m.  Substituting  these 
expansions  into  the  transformed  differential  equations  and  equating 
terms  of  like  order  of  ?  results  in  a  sequence  of  ordinary  differential 
equations.  In  Reference  5  the  solution  is  carried  to  first  order  (linear 
in  O  and  the  resulting  equations  are  given  here.  The  zero-order  equations 
are: 


f '  "  +  (n  -  ())  f  )f "  =  0 
o  o  o  o 


Pr 


+  =  -  y 


with  boundary  conditions: 


(B-3) 


f^(0)  =  f’(0)  =  0,  g^(0)  =  0 

f^(«)  =1,  g^(“)  =  1 


(B-4) 


The  prime  denotes  differentiation  with,  respect  to  n.  All  other  variables 
are  defined  in  other  sections  of  the  report.  The  first-order  equations, 
which  are  coupled  with  the  zero-order  equations,  are: 


£■"  +  (n  -  +  2«y„)£”  +  2(*^£;  -  i)£;  - 


a  F  /  o 
o ' 

<(>1 


o  o 


3<(>  f..  -  I  —  +  ijn 


f ' 
o 


/  <P1  \  ■ 

*  ^  r ) -  "♦d'o*  ■■  rt- 

'O'  o  o 

^  gl'  +  (n  -  ♦„£„)gi  £■  -  Dg^ 


Pr  I  a  F  ' 
\  o 


o  o  o  1 


a 


+  l)n  I  g’ 


(B-5) 


+  2 


Y  ^‘*’0^0  -^^Vr  yF  /  Y  («  <J>  2R  f'j^^o 

L  O  O  J  o  '  \  o  o  / 
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with  boundary  conditions: 


f^(0)  =  =  f’(»)  =  0 

g^(0)  =  g^(“>)  =  0 


(B-6) 


Equations  (B-3)  to  (B-6)  have  been  solved  for  and  Pr=0.72  and 

for  various  values  of  a  and  a.  The  resulting  solutions  for  f  and  f' 

o  1 

are  presented  in  Figure  B-1  and  the  solutions  for  and  g^  are  given  in 
Figure  B-2. 

Because  of  the  nature  of  the  solution  technique,  this  first-order 
solution  can  be  strictly  applicable  only  very  near  the  shock.  The  actual 
region  in  which  it  is  accurate  is  questionable  although  some  measure  can 
be  obtained  by  examining  the  accuracy  of  the  two  term  expansion  of  the 
external  flow  properties.  The  results  of  Appendix  A  indicate  that  the 
flow  properties,  in  particular  density,  change  very  rapidly  near  the 
shock.  The  two  term  expansion  of  the  density  function  R(C)  is  given  by 
extending  a  straight  line  from  the  shock  front  (5=0)  with  a  slope  equal 
to  the  slope  of  R  at  5=0.  By  this  straight  line  approximation,  R=0  for 
5=(y^“1)/[3(y+1)  +  o(y+5)].  This  yields  5=0.048  for  y=1»4  and  a=2.  The 
inviscid  flow  density  is  obviously  in  serious  error  for  5  less  than  that 
and  the  boundary  layer  solution  cannot  be  expected  to  be  any  more  accurate. 
Thus  the  boundary  layer  solution  presented  in  this  appendix  is  probably 
accurate  for  5<0.02  at  best. 
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APPENDIX  C 

METHOD  FOR  TURBULENT  BOUNDARY  LAYER  SOLUTION 


A  method  of  solving  the  transient  two-dimensional  equations  that 
govern  the  turbulent  boundary  layer  flows  generated  by  plane,  cylindrical, 

g 

and  spherical  strong  blast  waves  has  been  formulated  by  Quan.  Since  a 
predominant  portion  of  the  blast  wave  boundary  layer  is  turbulent,  this 
analysis  constitutes  a  logical  complement  to  the  laminar  boundary  layer 
study.  In  the  turbulent  flow  study,  the  integral  momentum  and  integral 
energy  equations  are  utilized  to  reduce  the  number  of  independent  variables 
from  three  to  two  (time  and  distance  along  the  surface).  A  similarity 
condition  further  reduces  the  two  independent  variables  to  one.  The 
resulting  ordinary  differential  equations  can  then  be  solved  by  standard 
numerical  procedure.  Quan's  formulation  of  the  solution  is  outlined 
in  this  appendix.  The  nomenclature  here  is  at  places  different  from 
that  of  the  main  text  of  this  report. 

The  transient  two-dimensional  compressible  turbulent  boundary  layer 

12 

equations  have  been  given  oy  Van  Driest.  It  can  be  shown  that  the 
equations  have  the  following  form  for  the  blast  wave  problem  where  the 
free-stream  flow  is  nonuniform  and  time-dependent: 


+  _L  .  9(pv)  Q 

9t  0  3x  9y 

X 


(C-1) 


+  -(pv)'h'  -  (pv)'u'  1^  (C-3) 
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where  the  bars  indicate  temporal  mean  values  and  the  primes  indicate 
instantaneous  fluctuations  from  the  mean.  The  symbols  p,  u,  v,  h,  and  y 
refer  to  the  density,  x  velocity  component,  y  velocity  component,  enthalpy, 
and  viscosity,  respectively;  p  is  the  pressure  at  the  edge  of  the  boundary 
layer;  Pr  is  the  Prandtl  number  which  is  assumed  to  be  constant;  and  t,  x, 
and  y  refer  to  time,  distance  along  the  surface,  and  distance  normal  to 
the  surface,  respectively.  The  temperature  is  given  by  T=h/c  where  c 
is  a  constant  specific  heat.  The  value  for  a  is  0  for  two-dimensional 
boundary  layer  and  1  for  axisymmetric  boundary  layer. 

There  are  two  common  methods  of  solving  turbulent  boundary  layer 
equations.  The  first  is  to  employ  seml-e.mpirical  expressions  relating 
the  fluctuation,  or  eddy  diffusivity,  terms  to  local  mean  properties, 
and  the  partial  differential  equations  governing  the  mean  values  are 
then  solved  by  finite-difference  procedure.  This  method,  at  least  in 
the  mathematical  sense,  is  quite  accurate.  However,  the  accuracy  and 
reliability  cf  the  semi-empirical  expressions  have  not  been  established 
at  present.  Moreover,  this  finite-difference  procedure  is  very  tedious 
even  for  steady-state  problems.  For  the  transient  problem  at  hand,  this 
method  is  unappealing,  to  say  the  least. 

The  second  method,  which  is  simpler  both  in  treating  the  physics 
of  turbulence  and  in  obtaining  numerical  results,  is  the  explicit  integral 
method.  Here,  the  differential  equations  are  integrated  across  the  boundary 
layer  and  the  dependency  on  the  y-direction  is  eliminated.  With  this 
method,  one  is  mainly  interested  in  the  gross  effect  of  turbulence  on  the 
boundary  layer  growth  and  not  in  the  detail  mechanics  of  turbulence. 

Instead  of  evaluating  the  local  diffusivity  terms  for  momentum  and  energy 
as  is  done  in  the  first  method,  one  prescribes  the  velocity  and  tempera¬ 
ture  profiles,  a  skin-friction  law  explicitly  relating  si.ear  stress  to 
momentum  thickness  and  an  explicit  relation  between  heat  and  momentum 
transfer.  The  resulting  partial  differential  equations  are  first-order 
and  they  describe  the  momentum  and  energy  thicknesses  as  functions  of 
X  and  t.  Two  questions  naturally  arise  with  the  integral  method  here. 

First,  are  the  boundary  layer  growth  results  sensitive  to  the  velocity 
and  temperature  shapes  that  are  assumed?  Second,  can  a  transformation 
be  found  such  that  x  and  t  are  related  by  a  single  variable? 
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Evidence  shows  that  the  answer  to  the  first  question  is  negative. 

In  fact,  the  integral  procedure  is  currently  the  standard  method  used  in 

14  15 

computing  turbulent  boundary  layers  in  nozzles  ’  where  finite  pressure 

gradients  exist.  It  is  also  used  in  shock  tube  turbulent  boundary  layer 
6  16 

studies  ’  in  which  pressure  gradients  are  absent.  It  should  be  noted, 
however,  that,  especially  for  the  present  problem,  the  skin-friction 
law  and  the  heat  and  momentum  correlation,  like  the  eddy  diffusivity 
expressions,  are  of  uncertain  validity. 

The  answer  to  the  second  question  is,  fortunately,  yes.  It  is  shown 
below  that,  indeed,  under  certain  restrictions  on  the  skin-friction  law, 
a  similarity  transformation  can  be  found.  Thus,  the  present  problem 
starting  with  second-order  differential  equations  with  three  independent 
variables  (t,  x,  and  y)  is  reduced  to  the  relatively  simple  task  of  solv¬ 
ing  first-order  ordinary  differential  equations.  To  our  knowledge,  a 
rigorous  transient  solution  to  two-dimensional  turbulent  flow  has  not 
been  obtained  previously  (of  course,  one  can  always  construct  what  may 
be  called  a  locally-similar  and  temporally-steady  solution  of  unassessed 
accuracy) .  The  method  of  combining  an  integral  technique  and  a  similarity 
transformation  appears  quite  attractive.  The  governing  equations  will  be 
derived  below  and  numerical  results  will  be  obtained  at  a  later  time. 

Equation  (C-1)  can  be  used  to  eliminate  pv;  Equation  (C-2)  can  be 
integrated  across  the  boundary  layer  to  yield: 


9 

T  =-;;—(puui)+p  . 

w  9t  e  e  e  9 


1  9  /  a  2  \ 

T7  **  +  ^  ^  (*  7 


9u 


+  p  u 


e  e  dx 


6*  (C-4) 


and  Equation  (C-3)  can  be  integrated  to  yield: 

9H 


q  =  |-  (p  H  f2)  +  p  + 

^w  3t  e  e  e  9t 


1  W  ^  u 
—  —  I  X  p  u  H 
a  9x  \  e  e  € 

X  \ 


9H 


(|)  +  p  u 


e  e  9x 


<S*  (C-5) 


where  x  is  the  shear  stress  at  the  wall  and  q  is  the  convective  heat 
w  w 

flux  to  the  wall.  The  subscript  e  denotes  conditions  at  the  edge  of  the 

2 

boundary  layer;  H  is  the  stagnation  enthalpy  (H=h+u  /2)  and 
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(C-6a) 

(C-6b) 

(C-6c) 

(C-6d) 

(C-6e) 

(C-6f) 


The  dispalcetnent  thickness  6*,  momentum  thickness  6,  and  energy  thick¬ 
ness  iji  are  those  commonly  encountered  in  steady-state  problems.  It  is  seen 
that  three  other  th''ckness  parameters  are  introduced  in  the  transient 
problem.  We  shall  temporarily  call  X*  the  density  thickness,  w  the 
acceleration  thickness,  and  0,  the  enthalpy  thickness. 

In  deriving  Equations  (C--4)  and  (C-5)  ,  it  has  been  assumed  that 

*00  ^oo 

pu  dy  ~  p  u  dy  (C-7a) 

o  •'  o 


(C-7b) 

(C-7c) 
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This  assumption  is  generally  employed  in  integral  methods  for  turbulent 
17 

flow.  Also,  the  following  boundary  conditions  have  been  employed; 


u=0,  v*0,  h=0 

at  y  =  0 

(C-8a) 

u  =  u  ,  h  =  h 
e  e 

at  y  ->■  w 

(C-8b) 

It  is  now  assumed  that  the  velocity  and  stagnation  enthalpy  can  be 
represented  by  power  law  profiles  of  the  form 


y  _<  6 

(C-9a) 

1 

u 

e 

y  >  6 

(C-9b) 

y  £  A 

(C-lOa) 

JL=  1 

H 

y  >  A 

(C-lOb) 

where  6  is  the  velocity  thickness,  A  is  the  temperature  thickness,  and  n 
is  an  empirical  number  usually  given  the  value  of  1/7. 


From  Equations  (C-9)  and  (C-10)  and  the  definition  of  H,  one  can 

pressl 
Integration 


obtain  p/p^  as  a  function  of  (y/6)  and  (y/A)  .  Its  explicit  expression 


depends  on  whether  y  A  or  y  >  A  and  whether  6  ^  A  or  6  >  A 
of  Equation  (C-6)  yields 

6*/6 

e/6 


f^(C) 


f2(?) 


X*/6 
0)/ 6 
U/6 


f^(?) 

f5(0 

fgU) 


(C-lla) 

(C-llb) 

(C-llc) 

(C-lld) 

(C-lle) 

(C-llf) 
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where  ?  =  A/(S.  Using  Equation  (C-11) ,  one  may  write  Equations  (C-4)  and 
(C-5)  in  the  form 


w 


8u 


+  p  u 


8,(0  0 


^w 


“  ^  Pe  if  ^  ^  (^'"pe^v) 


e  e  8x  “1 
(C-12) 
3H 


4-  p  u 


e  e  9x  ®4 
(C-13) 


where 

gfC^) 

82(5) 

83(4) 

= 

83(4) 
8^(4)  = 


f^(c)/f2(0 

(C-14a) 

f4(4)/f2(4) 

(C-14b) 

f5(C)/f2(C) 

(C-I4c) 

f3^(c)/f3(4) 

(C-I4d) 

f4(C)/f3(c) 

(C-14e) 

(C-14f) 

It  is  noted  that  A  can  be  written  as  a  function  of  (j)/e  through 
Equations  (C-llb)  and  (C-llc) .  Equations  (C-12)  and  (C-13)  can  then  be 
solved  for  6  and  if>  simultaneously  provided  that  relations  for  and  q^ 
as  functions  of  6  and  <fi  are  given.  We  are,  however,  interested  in  a  more 
ambitious  step,  namely,  in  reducing  Equations  (C-12)  and  (C-13)  to 
ordinary  differential  equations.  In  this  regard,  we  note  that  similarity 
exists  in  the  free-stream  conditions  for  strong  shocks  for  which 


where 


Poo^^(5) 

(C-15a) 

Ug(T))i)(5) 

(C-15b) 

uJ(T) 

Y  R({)  +  2 

H  uJ(T)e(u 

(C-15c) 

C  =  1  -  — 
ct“ 

(C-16a) 

T  =  t 

(C-16b) 
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and  p  ,  Y>  C,  and  m  are  constants,  and 


u  (t)  =  Cmi 
s 


m-l 


(C-17) 


Equation  (C-12)  can  be  written  as 


+  p^Rg2X  6 


Ba-ou 

T  a  d5  dT 


m  3 5 


2  3 


Ct 


CC-18) 


We  now  assume  that  ?  or  the  ratio  iji/6  depends  only  on  ^  and  that 
6  can  be  written  in  the  form 


e(x,t) 


(C-19a) 


or 


x'^eCx.t)  =  BT®e^(5) 


(C-19b) 


where  the  constant  a  is  to  be  determined  and  the  constant  B  can  be 
arbitrarily  selected  for  convenience.  Equation  (C-18)  becomes 


1  a  m  3d 

— ;r  T  X  ==  —  (l-5)u  T  — 

p  3  w  T  s  dC 


+  Rg2T  6^ 


,  .  du  1 

—  (l-C)u  -jy  +  'p-r— 
T  s  d5  di 


1  2  a  d  4,,2„  \  1.2,  a„  d)/; 

- u  T  -rr  (R'i'  0,  ) - Ru  <pg,T  0, 


(C-20) 
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Equation  (C-20)  can  be  written  in  the  form 


nH-a-2-m0„ .  _ . 
T  =  T  GiO 

w 


(C-21) 


where 


p^BC^  ‘^m  ( 

G(0  = 


(1-0 


+  (m-l+a)Ri|jg20j^ 


+  (m-l)iJjJ 


dij; 

dK 


(C-22) 


The  question  now  is  can  one  find  an  expression  for  that  is  a 
'ur.ctior;  ci  and  varies  with  t  to  some  power.  Fortunately,  a  physically 
a:ce?ta;I«.  expression  tor  such  purpose  does  exist.  For  turbulent  flow 
in  noz;les‘°  and  In  shock  tubes, a  relation  between  skin-friction 
and  hourcary  layer  thickness  can  be  taken  as 


where 


and 


(C-23) 


(C-24) 


T  =  0.5  (T  +T  )  +  0.22  (T  -T  )  (G-25) 

m  we  re 

Assuming  y  to  vary  directly  with  T,  considering  T  <<  T  ,  taking  the 

1/3  2  ^  ^ 

recovery  temperature  to  be  using  h=CpT,  one 

obtains 
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and 


i  = 


-,1/2 


1/3  2 

0  5  ^  O.llYPr 


j(5) 


(C-26) 


^  ^  ^2m-2+(in-l+ma-a)/4 


(C-27) 


where 


2  9/4 

K(C)  =  0.0225  j  pM 


YFc'^(l-5)'^y„f2 


(Y-l)B,i;p„Rh^6^ 


1/4 


(C-28) 


Equating  the  powers  of  t  in  Equations  (C-21)  and  (C-27) ,  one  obtains 


m+a-'2-ma  =  2m-2+(in-l+ma-a) /4 


(C-29) 


which  yields 


a  “  m(a+l)  -  1/5 


(C-30) 


Thus,  a  similarity  solution  can  be  achieved  for  0.  It  may  be  pointed 
out  here  that,  while  other  power-law  expressions  for  still  allow  for 
similarity  transformation,  certain  nonpower-law  expressions,  which  may  be 
more  accurate,  do  not  yield  similarity  conditions.  The  sacrifice  of  such 
possible  gain  in  accuracy  is  of  r  .  concern  here  since  turbulent  flow  is 
subject  to  a  great  many  other  inaccuracies  and  one  is  interested  in  obtain¬ 
ing  a  reasonable  solution. 


For  convective  heat  transfer.  Equation  (C-13)  can  be  written  as 

du 


m  ,  .  2  dg  . 

—  (l-C)u  -jT  +  23u  -7 

T  '  s  d5  s  d 


i 

T 


+  't>| 

- —  p  u  ^  - p  Ru  (C-31) 
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Let  (ti(x,t)  be  given  by 


x'^4i(x,t)  = 


Then  Equation  (C-31)  can  be  written  in  the  form 


2tiH-b-3-ma, 

%  =  -r  L(5) 


(C-32) 


(C-33) 


where 

p  BC^  “^m^  (  ,  /  \ 

L(?)  =  Jtn(l-5)  +  (2in-2+a)R3gg^^ 

+  j^m(l-C)  II  +  2(m-l)3 

-  m  ^  -  mRiJ^g^4>i  ||  |  (C-34) 

It  is  now  needed  to  find  an  expression  for  that  varies  with  t  to  some 
power.  In  addition,  this  expression  must  yield  a  value  of  b  that  is 
equal  to  a  in  order  for  the  assumption  of  <^/Q  as  a  function  of  only  C  to 
hold.  For  this  purpose,  one  may  employ  an  expression  that  has  been  used 
in  shock  tube  studies 

c  T  _ 

q  =  -2-^  Pr"^'^  T  (C-35) 

^w  u  w 

e 

which  yields 


where 


3m-3+(m-l-hna-a)  74^ 

T  Q(5) 


(C-36) 


Q(5)  =  CmPr 


-2/3 


/izl  JL  +  1 

I  Y  m  2 


l/3„ 


(C-37) 
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Equating  the  powers  of  t  in  Equations  (C-33)  and  (C-36)  and  using 
Equation  (C-30)  ,  it  is  found  that  indeed 

b  =  a  =  m(a+l)  -  1/5  (C-38) 

Thus,  a  similarity  solution  can  also  be  obtained  for  (j).  Equating 
G(5)  to  K(5)  and  equating  L(5)  to  Q(C)  results  in  two  ordinary  differential 
equations  to  solve  for  0^  and  The  calculations  will  be  carried  out 

in  detail  at  a  later  date. 

In  summary,  the  equations  that  govern  blast  wave  turbulent  boundary 
layer  flows  are  second-order  partial  differential  equations  with  three 
independent  variables.  In  the  present  report,  the  integral  method  is 
used  to  reduce  the  equations  to  first-order  and  the  number  of  independent 
variables  to  two.  Similarity  transformation  is  shown  to  be  achievable 
and  is  used  to  reduce  the  problem  to  the  relatively  simple  task  of  solving 
two  coupled  ordinary  differential  equations.  These  equations  are  derived 
in  the  present  report  and  are  to  be  progreiomed  for  numerical  results. 
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ABSTRACT 

This  second  part  of  the  report  is  on  particle  entrainment.  Using  a  lift 
force  generated  by  the  boundary  layer  velocity  gradient,  the  particle  velo¬ 
city  is  estimated.  Sample  particle  trajectories  in  strong  blast  wave  flow- 
fields  are  illustrated,  and  the  amount  of  soil  erosion  due  to  aerodynamic 
entrainment  is  assessed. 
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NOMENCLATURE 

a  =  particle  radius 

A  =  particle  cross-sectioniil  area;  soil  surface  area 
Cp  =  specific  heat  of  gas  at  constant  pressure 
c  =  specific  heat  of  particle 

C  =  constant  related  to  explosion  strength 
=  drag  coefficient 
C  =  lift  coefficient 

d  =  particle  diameter 

D  =  particle  diameter;  pipe  diameter 

f  =  dimensionless  average  momentum  of  entrained  particles 

F  =  function  given  by  Equation  (3-9) 

Fj^  =  lift  force 

F  “  horizontal  drag  force 

X  ° 

g  =  gravitational  acceleration 

G  =  function  given  by  Equation  (3-10) 

k  =  gas  conductivity 

K  =  3p/8pg 

i  =  boundary  layer  parameter 

=  boundary  layer  parameter 
L  =  heat  of  vaporization  for  particle 

m  =  mass  of  a  particle;  mass  of  soil  per  unit  surface  area 
M  =  total  mass  of  erosion 

Nu  =  Nusselt  number  based  on  particle  diameter 
p  ==  pressure  in  gas 

=  pressure  in  soil 
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NOMENCLATURE  (Continued) 

ag/K'-i  ^ 

ti 

Reynolds  number  based  on  particle  diameter 
C^6/a 

time 

gas  temperature 

particle  temperature 

vaporization  or  ablation  temperature 

gas  velocity  in  x-direction 

particle  velocity  in  x-direction 

gas  velocity  in  y-dlrection 

particle  velocity  in  y-direction 

relative  velocity  between  gas  and  particle 

distance  along  soil  surface 

distance  normal  to  soil  surface 

boundary  layer  thickness 

height  reached  by  a  particle  due  to  direct  lift 
a/ 6 

coefficient  for  soil  erosion 

gas  specific  heat  ratio 

boundary  layer  parameter,  =  y^/4.6 

boundary  layer  parameter 

gas  viscosity 

gas  dynamic  viscosity 

gas  density 

density  of  a  particle 
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NOMKNCLATURK  (Continued) 

T  =  shear  stress  on  soil  surface 

o 

T*  =  restraining  stress  of  soil 

Subscripts 

e  =  at  edge  of  boundary  layer  during  ascent 

f  =  at  edge  of  boundary  layer  during  descent 

g  =  at  ground  position 

i  =  initial  condition 

m  =  at  maximum  height  of  a  particle 

p  =  particle 

s  =  particle 

w  ==  wall  condition 

00  =  ambient  condition 
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I.  INTRODUCTION 

In  a  nuclear  explosion,  soil  can  be  lifted  from  the  ground  by  several 
mechanisms  such  as  crater  splash,  vaporisation,  elastic  rebound,  etc. 

These  processes  are  inter-dependent  and  a  complete  study  of  the  simultaneous 
processes  has  not  been  performed.  Not  only  is  the  physics  of  the  lofting 
mechanisms  not  yet  well  understood,  but  the  rigorous  mathematical  models 
that  can  be  proposed  are  very  difficult  to  solve.  The  present  study  is 
only  concerned  with  the  aerodynamic  effects  on  particle  entrainment,  and 
the  interaction  effects  of  other  mechanisms  are  not  considered. 

On  the  soil  erosion  problem,  one  is  mainly  interested  in  answers  to 
two  questions:  what  is  the  initial  velocity  of  a  particle  as  it  leaves  the 
ground  (or  as  it  leaves  the  boundary  layer)  and  how  much  dust  of  various 
particle  sizes  is  entrained?  Answers  to  these  questions  would  allow  one  to 
predict  the  dust  density  distribution  as  a  function  of  space  and  time  for 
a  given  gas  flowfield. 

The  aerodynamic  ertrainment  of  surface  particles  into  free  streams 

has  been  investigated  analytically  and  experimentally  by  workers  of  various 

1  2  3 

disciplines.  Early  studies,  *  as  well  as  some  recent  ones,  were  con¬ 
cerned  with  wind  erosion  of  soil.  About  a  decade  ago,  the  aerodynamic 
forces  on  surface  particles  were  Investigated  in  connection  with  the  VTOL 
aircraft  and  helicopters  downwash  problem.^  About  the  same  time,  studies 
were  made  to  assess  the  dust  entrainment  caused  by  the  rocket  exhaust  of  a 
spacecraft  during  lunar  landing.^ 

The  investigations  indicated  above  represent  the  pioneering  efforts  in 
their  respective  fields,  and  as  such  they  constitute  significant  contribu¬ 
tions.  However,  their  results  do  not  appear  to  be  directly  applicable  to 
or  sufficient  for  the  problem  of  predicting  the  dust  entrainment  caused  by 
a  blast  wave.  The  wind  erosion  studies  are  primarily  one-dimensional,  i.e., 
the  wind  variation  along  the  surface  and  with  time  is  neglected  while  a 
blast  wave  flowfield  is  strongly  position-  and  time-dependent.  The  downwash 
impingement  results  of  Vidal^  provide  an  estimate  of  the  lift  and  drag 
forces  of  a  particle  when  it  is  on  the  ground,  but  the  subsequent  motion 
of  the  particle  is  not  considered.  For  dust  entrainment  due  to  lunar 
landing,  Roberts^  assumed  that  the  aerodynamic  shear  stress  on  the  surface, 
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minus  some  value  which  represents  a  restraining  stress,  is  proportional 
to  the  rate  of  transfer  of  momentum,  per  unit  area,  to  the  particles.  The 
proportionality  factor  is  related  to  the  particle  size  and  rocket  exhaust 
conditions  and  is  not  directly  applicable  to  blast  wave  conditions.  Some 
of  the  results  and  ideas  of  Vidal  and  Roberts  will  be  employed,  however,  in 
the  present  study. 

For  the  problem  of  aerodynamic  entrainment  of  soil  by  a  blast  wave, 
several  studies  have  been  performed  recently.  Swatosh  and  Wiedermann^ 
assumed  both  the  mass  erosion  rate  and  the  vertical  velocity  (considered 
to  be  due  to  turbulence  in  the  flow)  to  be  proportional  to  the  free- 
steam  air  velocity;  and  they  performed  experiments  to  determine  the  propor¬ 
tionality  factors.  Their  results  are  important  contributions  to  the  present 
knowledge.  However,  since  these  results  are  not  directly  related  to 
boundary  layer  properties,  their  applications  are  probably  limited  if  shear 
stress  is  an  important  factor  in  determining  erosion.  For  example,  consider 
uniform  flow  over  a  flat  plate.  Since  the  shear  stress  on  the  surface 
decreases-  rapidly  with  distance  from  the  leading  edge,  it  is  reasonable 
to  expect  a  higher  erosion  rate  near  the  leading  edge.  However,  the  rela¬ 
tion  given  by  Swatosh  and  Wiedermann  would  show  a  uniform  erosion  rate. 

The  fact  that  boundary  layer  properties  probably  play  an  important  role  in 
determining  soil  erosion  may  partially  account  for  the  variations  in  the 
experimental  data.  For  example,  the  erosion  constant  for  variable  air 
velocity  is  generally  much  higher  (sometimes  by  two  orders  of  magnitude) 
than  that  for  steady  air  velocity;  this  may  be  due  to  effects  related  to 
boundary  layer  build-up, 

Trulio  and  others^  employed  a  simple  model  by  assuming  that  the  hori¬ 
zontal  impulse  delivered  to  the  ground  (due  to  shear  stress)  is  entirely 
converted  into  vertical  momentum  of  particles.  The  particles  are  considered 
to  rotate,  slide,  or  bounce  along  the  surface  until  they  accelerate  by  the 
horizontal  aerodynamic  drag  force  to  a  sufficient  speed  to  bounce,  or  dis¬ 
lodge  other  particles,  off  the  surface  and  into  the  free  stream.  The 
particles  are  considered  to  leave  the  ground  vertically  with  a  speed 
inversely  proportional  to  the  particle  diameter;  the  maximum  speed  is  taken 
to  be  half  the  free-stream  gas  speed.  Due  to  the  lack  of  a  rigorous 
boundary  layer  solution,  a  local  one-dimensional  steady  turbulent  boundary 
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layer  profile,  which  has  been  verified  to  be  approximately  correct  for 
sandstorms,  is  employed  to  evaluate  the  shear  stress.  Again,  the  blast 
wave  flowfield  is  actually  highly  transient  and  contains  severe  spacial 
gradients  in  properties.  The  validity  of  using  any  steady  and  one¬ 
dimensional  model  is  questionable.  It  is  the  desire  of  obtaining  a  rigorous 
boundary  layer  solution  that  led  to  the  current  boundary  layer  studies  by 
TRW  Systems. 

In  the  present  study,  the  particle  motion  near  the  surface  due  to 
aerodynamic  forces  is  estimated.  Sample  particle  trajectory  calculations 
are  made.  The  amount  of  soil  entrained  is  assessed  using  the  laminar 
boundary  layer  solution  developed.  Since  the  blast  wave  boundary  layer  is 
expected  to  be  predominantly  turbulent,  a  turbulent  flow  solution  is  being 

O 

pursued.  Wlien  this  solution  is  completed,  the  amount  of  soil  entrained  can 
be  similarly  assessed. 
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2.  PARTICLE  LIFT  IN  BOUNDARY  LAYER 


The  existence  of  a  velocity  difference  between  a  particle  and  a  fluid 

results  in  an  aerodynamic  force  between  them.  The  force  component  in  the 

direction  of  the  relative  velocity  is  known  as  the  drag  and  the  component 

normal  to  it  as  the  lift.  The  drag  for  a  spherical  particle  in  uniform, 

steady  incompressible  flow  is  a  function  of  only  the  Reynolds  number  and 

9 

is  well  known  up  to  a  Reynolds  number  of  about  a  million.  The  effects  of 
particle  shape,  flow  unsteadiness,  and  compressibility,  etc.  on  drag  have 
been  investigated  (e.g..  References  10,  11,  and  12),  and  some  semi-empirical 
correlations  are  available.  However,  the  lift  on  a  sphere  produced  by  a 
nonuniform  flowfield  such  as  a  boundary  layer  is  still  not  well  understood, 
although  it  has  been  postulated  to  be  the  responsible  cause  for  particle 
migration  from  a  lower  to  a  higher  velocity  region.  An  indication  of  the 
complexity  of  the  lift  mechanisms  is  the  fact  that  the  force  on  a  rotating 
sphere  in  a  uniform  fluid,  which  is  due  to  what  is  known  as  the  Magnus 
effect  and  which  accounts  for  the  irregular  flight  of  a  tennis  ball  or 
baseball,  has  been  studied  for  about  three  centuries  and  the  magnitude  of 
this  force  is  still  undetermined. 

Recently,  there  appeared  several  analyses  dealing  with  the  lift 

13 

exerted  on  a  spherical  particle  by  a  shear  flow.  Eichhorn  and  Small  per¬ 
formed  experiments  for  spheres  suspended  in  Poiseuille  flow  and  obtained 
the  following  relation 


=  7  X  10 


(2-1) 


where  C  is  the  lift  coefficient  (the  ratio  of  the  lift  force  to  the  product 
of  the  sphere's  cross-sectional  area  and  the  fluid  dynamic  pressure),  and 
d,  u,  y,  D,  and  Re  refer  to  particle  diameter,  fluid  velocity  at  center  of 
particle,  distance  from  pipe  wall,  pipe  diameter,  and  particle  Reynolds 
number  respectively.  The  experiments  were  performed  using  a  0.419  in.  dia¬ 
meter  tube  and  for  d  ranging  from  0,061  in  to  0.126  in..  Re  from  80  to  250, 
and  (d/u)(du/dy)  from  0  to  1.1.  The  value  of  C  obtained  range  from  0  to 
about  1.0. 
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While  Equation  (2-1)  showed  C  to  vary  with  the  velocity  gradient 

14 

(du/dy)  to  the  2.0  power,  Saffman  found  that  at  low  particle  Reynolds 
number  C  varies  with  the  gradient  to  0.5  power: 

Li 

where  K  =  6. 46,  v  is  the  fluid  kinematic  viscosity,  and  the  distance  y  is 
measured  normal  to  the  flow  direction.  Saffman* s  analysis  is  for  uniform 
shear  flow  and  restricted  to  particle  Reynolds  numbers  less  than  1.0. 

4  15 

Vidal  ,  by  using  Hall's  results  for  the  tangential  velocity  variation 
on  a  sphere  in  uniform  shear  flow,  found  that 

=  0.998  k  (2-3) 

where 

k  .  ^  (2-4) 

u  dy 

and  a  is  tne  particle  radius.  Thus,  Vidal's  relation  shows  to  vary 
linearly  with  (du/dy).  Equation  (2-3)  becomes  Invalid  for  k  1. 

A  purpose  of  the  present  study  is  to  illustrate  the  effects  of  lift  on 
dust  entrainment.  Since  there  does  not  exist  an  established  general 
formula  to  predict  lift,  it  is  advantageous  to  choose  a  least  restrictive 
relation  for  the  purpose  of  illustration.  For  a  particle  in  a  high  speed 
boundary  layer,  the  Reynolds  number  is  very  large  (compared  to  1.0)  except 
for  small  particles  at  or  near  the  ground.  Small  particles,  however,  have 
large  values  of  drag  coefficient;  even  if  they  are  lifted  off  the  surface, 
they  soon  approach  and  follow  the  fluid  streamlines  and  will  not  reach  a 
significant  vertical  distance  above  the  ground.  Therefore,  Saffman's 
results  will  not  be  utilized  in  the  present  study.  The  pipe  flow  experi¬ 
ment,  though  it  can  conceivably  be  applied  to  external  flows  by  relating  , 
the  pipe  diameter  to  boundary  layer  thickness,  is  also  limited  in  particle 
Reynolds  number  range  and  in  the  ratio  of  particle  size  to  boundary  layer 
thickness.  Therefore,  the  relation  for  lift  coefficient  given  by  Eichhorn 
and  Small  will  also  not  be  employed.  Vidal's  relation.  Equation  (2-3),  is 
used  in  the  present  study.  This  relation  is  simple  to  apply,  and  its 
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restriction  of  k  <1,  as  will  be  shown  later  in  this  section,  is  not  a 
severe  limitation.  It  should  be  noted,  however,  that  Vidal's  relation  is 
based  on  uniform  shear  flow  over  a  particle  at  rest.  In  the  present  study 
of  particle  motion  in  a  boundary  layer,  this  relation  is  applied  by  using 
a  local  velocity  gradient  and  a  local  relative  velocity. 

The  particle  trajectories  inside  and  outside  the  boundary  layer  in 
the  transient  nonuniform  blast  flowfields  will  be  shown  in  the  next  section. 
In  the  present  section,  it  is  desired  to  acquire  some  physical  insight  and 
some  qualitative  measure  of  the  effects  of  various  parameters  on  the 
potential  motion  of  a  particle  as  it  is  lifted  off  a  surface.  For  this 
purpose,  a  steady  one-dimensional  boundary  layer  (i.e.,  a  boundary  layer 
of  constant  thickness  along  a  surface)  with  constant  properties  is 
considered.  The  velocity  distribution  is  represented  by 

u  =  u  (1  -  (2-5) 

e 

V  =  0  (2-6) 

where  u^,  u,  v,  and  y  denote  the  gas  free-stream  velocity,  velocity  parallel 
to  the  ground,  velocity  normal  to  the  ground,  and  distance  from  the  surface, 
respectively.  The  parameter  6  is  a  measure  of  the  boundary  layer  thickness; 
at  y  =  4.6  ,  the  velocity  becomes  u  =  0.99  u^. 

The  particle  is  considered  to  be  spherical,  and  u  and  y  are  measured 
at  the  particle  center.  It  is  assumed  that  the  particle  is  initially  at 
rest  and  tangent  to  the  ground,  i.e.,  y  =  a  initially.  The  equations 
governing  the  particle  motion  are  taken  as 

'^p  ^  =  2  K  ■  ^p^^  +  Cp  V  (V  -  Vp)]  -  mg  (2-7) 

du 

m  V  =  -;rPAC,^V(u-u)  (2-8) 

p  dy  2  D  p 

where  m.  A,  p,  C^,  g,  u^ ,  v^,  and  V  denote,  respectively,  the  particle  mass, 
particle  cross-sectional  area,  fluid  density,  drag  coefficient,  gravitational 
acceleration,  particle  velocity  parallel  to  the  ground,  partical  velocity 
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normal  to  the  ground,  and  relative  velocity  between  particle  and  fluid.  In 
deriving  Equations  (2-7)  and  (2-8),  it  has  been  assumed  that  lift  is  impor¬ 
tant  only  in  the  vertical  direction  and  that  the  lift  coefficient  is  based 
on  relative  horizontal  fluid  and  particle  velocities.  Considering  Equations 
(2-3)  and  (2-4),  one  may  take 


C 


L 


a  du 

u  -  u  dy 
P 


(2-9) 


where  the  coefficient  0,998  has  been  replaced  by  1.0  for  simplicity  and 
where  u  ^  du/dy  has  been  interpreted  as  (u-u  )  ^du/dy  for  a  particle 
in  motion.  Equations  (2-7)  and  (2-8)  with  v  =  0  can  be  written  as 


^  ,  du  „ 

(u-u)^ - Vv 

P  dy  a  p 


(2-10) 


where 


K  —  V  (u  -  u  ) 
a  p' 


K 


3P 

8p 

s 


(2-11) 


(2-12) 


V  =  [(u  -  Up)^  +  (2-13) 

and  where  P^  denotes  the  particle  material  density. 

At  this  point,  it  appears  advantageous  to  examine  the  variation  of 

the  lift  force  From  Equation  (2-10)  it  is  seen  that  Fj^  is  proportional 

to  (u-u  )du/dy.  Within  a  short  distance  from  the  surface,  u  is  negli- 
P  P 

gible  compared  to  u,  and  the  lift  force  is  proportional  to  u(du/dy)  which 

varies  with  y  according  to  Equation  (2-5)  so  that 


m 


(1  - 


(2-14) 


2 

Thus,  the  lift  force  per  unit  mass  has  a  maximum  value  of  Ku  /46  and 

2  ^ 

occurs  at  y  =  0.6936.  The  variation  of  F  6/mKu  with  y/6  is  shown  in 

1.4  6 

Figure  2-1. 
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From  Equation  (2-10),  it  is  seen  that  K(udu/dy)  at  y  =  a  must  be 
greater  than  g  in  order  for  a  particle  that  is  initially  resting  on  a  flat 
surface  to  be  lifted.  Using  Equation  (2-5),  this  requirement  can  be 
written  as 

(Ku^^/ag)C  >  1  (2-15) 

where 

5  =  (a/6)(l  -  (2-16) 


The  function  5  is  shown  in  Figure  2-2  and  has  a  maximum  value  of  0.260  at 
a/6  =  1.45.  The  interesting  point  is  that  for  a  given  particle  size,  there 
exists  a  preferential  range  of  boundary  layer  thickness  for  initial  entrain¬ 
ment.  That  is,  lift  may  not  exist  near  the  leading  edge  on  a  surface 
where  the  boundary  layer  is  thin  (although  the  shear  stress  is  very  high 
there)  or  near  the  tail  end  where  the  boundary  layer  is  thick.  If  the 
flow  conditions  are  such  that  Equation  (2-17)  is  not  satisfied  everywhere 
along  a  surface,  then  no  lift  may  occur  at  all.  Also,  in  order  for  the 

particle  to  be  lifted  to  a  significant  distance  above  the  ground,  the 
2 

condition  of  (Ku^  /ag) C  >>  1  is  required. 

It  is  interesting  to  compare  the  lift  force  Fj^  with  the  horizontal 
drag  force  F  .  From  Equations  (2-10)  and  (2-11) ,  the  ratio  is  given  by 


1  ^  ^ 
Fx  “  Cjj  V  dy 


(2-17) 


For  an  estimate  of  the  value  of  this  ratio,  one  may  approximate  V  by  u  to 
obtain 


X 


^  ^  -  k  -  _ 


1) 


(2-18) 


The  parameter  k  is  plotted  as  a  function  of  a/6  and  y/a  in  Figure  2-3.  It 
is  seen  that  k  has  a  maximum  limiting  value  of  1.0  which  occurs  when  the 
particle  is  at  the  surface  (y  =  a)  and  when  a/6  0.  Thus,  the  condition 

of  k  <  1.0  as  required  for  the  lift  relation  [Equations  (2-3)  and  (2-4)]  is 
satisfied.  (However,  if  the  particle  is  partially  imbedded  in  the  ground 
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initially,  i.e.,  y  <  a,  then  k  >  1.0  initially  except  for  large  values 
of  a/ (5.)  At  a  distance  of  a  few  radii  away  (say  y/a  =  10),  the  ratio  of 
vertical  lift  to  horizontal  drag  becomes  negligible. 


Attention  is  now  directed  towards  the  solution  of  Equations  (2-10)  and 
(2-11).  From  these  equations,  it  can  be  deduced  that  a  particle  is  first 
accelerated  off  the  ground  by  the  lift  force,  avid  then  decelerated  verti¬ 
cally  by  drag  and  gravity.  During  its  subsequent  return  towards  the  ground, 
the  particle  may  get  lifted  again  if  the  lift  at  any  point  within  the  bound¬ 
ary  layer  overcomes  gravity  and  a  few  such  damped  oscillations  may  occur. 
However,  during  its  descent,  a  particle  generally  has  attained  a  horizontal 
velocity  nearly  equal  to  the  fluid  velocity  and  the  lift  due  to  (u-u  ) (du/dy) 
is  negligible.  In  fact,  u^  becomes  greater  than  u  somewhere  within  the 
boundary  layer,  and  a  negative  lift  develops  which  aids  gravity  in  returning 
the  particle  to  the  ground. 


In  order  to  solve  Equations  (2-10)  and  (2-11)  analytically,  it  is 
necessary  to  make  some  simplifications.  The  drag  coefficient  will  be 
taken  as  a  constant;  this  is  acceptable  since  the  particle  Reynolds  number 
is  generally  very  large.  The  flow  is  divided  into  two  regions  for  solution: 
vertical  acceleration  of  the  particle  within  the  boundary  layer  and  vertical 
deceleration  outside  the  boundary  layer. 


The  boundary  layer  thickness  will  be  denoted  by  y^  and  has  a  value  of 

4.66  at  which  u  =  0.99  u^.  For  acceleration  in  the  boundary  layer,  it  is 

assumed  that  vertical  drag  and  gravity  are  negligible  compared  to  lift. 

Particles  that  are  too  small  to  satisfy  this  assumption  will  not  be  lifted 

significantly  above  the  boundary  layer  and  will  not  be  considered.  It  is 

also  assumed  that  u^  is  negligible  compared  to  u  in  computing  lift.  These 

assumptions  are  somewhat  conservative,  i.e.,  they  would  result  in  over- 

predicting  the  values  of  the  vertical  velocity.  The  boundary  conditions 

should  be  Vp  =  Up  =  0  at  y  =  a  if  a  particle  is  initially  at  rest  on  a  flat 

surface.  However,  the  boundary  will  be  given  an  arbitrary  value  of  y  =  y^ 

instead  of  y  =  a  in  the  solution.  With  these  assumptions  and  boundary 

conditions,  the  solution  for  v  is 

P 


V 

P 


(Ku^  -  Ku^^)^'^^ 


for  y  y^ 


(2-19) 
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where  u  is  given  by  Equation  (2-5)  and  is  u  at  y^.  Equation  (2-19) 
may  be  written  as 


-  =  2e  ^(l-e^)-e^^(l-e  ^^)  for  y  £  y  (2-20) 

Ku  ® 

e 


where 


^  =  y^/iS 


z  =  (y  -  y^)/<5 


(2-21) 

(2-22) 


Equation  (2-19)  shows  that  the  maximum  attainable  value  of  v  due  to  lift 

1/2  ^ 

is  only  a  small  fraction  of  u  ,  i.e.  ,  v  K  u  where  K  is  on  the  order 
_3  e’  ’  p  e 

of  10  (when  the  fluid  is  a  liquid,  however,  K  is  of  order  1).  Also, 

Equation  (2-19)  shows  that  v^  is  smaller  for  larger  particles  since  u^  is 

larger  for  larger  particles. 

In  using  Equation  (2-11)  to  solve  for  u  in  the  boundary  layer,  V  is 

1/2  ^ 

approximated  by  u  and  v^  by  K  u.  The  result  is 


S.  = 


u 


=  1  -  e 


-sz  .  s  -  X  ,  -z  -sz 


1-s 


e  (e 


0 


(2-23) 


where 


s  =  Cjj/a 


(2-2A) 


«  =  a/6 


(2-25) 


If  s  =  1,  the  solution  for  u^  can  be  obtained  by  applying  L'Hospital's 
rule  to  Equation  (2-23) . 


The  time  a  particle  takes  to  travel  from  y,  to  y  is  given  by 

/  y  -1  ^ 

=  Vp  dy  and  the  horizontal  distance  x  a  particle  has  traveled  is 


given  by  X  =  (u^/Vp)  dy. 
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After  the  particle  has  reached  the  edge  of  the  boundary  layer,  no 

lift  due  to  shear  exists  and  the  particle  is  considered  to  be  decelerated 

vertically  by  drag  and  gravity.  In  this  region,  V  can  be  approximated  by 

u  -u  .  It  is  convenient  to  choose  a  transformed  time  coordinate  4 
e  p 

instead  of  z  as  the  independent  variable  in  this  region.  The  solution  is 


V  ,  /v 
%  ^ 


2Cj^w 


for  t  >  t 

e 


(2-26) 


where  e  denotes  values  at  y  and 

■^e 


r 


(2-27) 


(2-28) 


w 


(2-29) 


41  =  1  +  KC_u  w  (t  -  t  )/a 
D  e  e 


(2-30) 


The  distances  y  and  x  are  related  to  t  by  y  =  y  +  j  ^  v  dt  and 

/t  ®  -7  to  p 

H  u  dt, 

tp  p  ’ 


X  “  X  +  /^  u  dt,  or 
e  y  te  p  ’ 


*  "  ^e  KCjjW 


for  t  >  t 
—  e 


(2-31) 


for  t  >  t 
—  e 


(2-32) 


The  maximum  height  occurs  at  v^  =  0  and  can  be  obtained  by  using 
Equations  (2-26)  and  (2-31) .  The  time  it  takes  a  particle  to  return  to 
the  ground  and  the  velocity  components  there  can  be  obtained  by  setting 
z  =  0  in  Equation  (2-31)  and  using  Equations  (2-26) ,  (2-27) ,  and  (2-30) . 
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Examination  of  Equations  (2-19)  to  (2-32)  shows  that  and  u  /u^ 

are  dependent  on  five  dimensionless  parameters:  K,  Cj^,  a  ,  r,  and  z  or  ((>  . 
(The  boundary  of  =  a  is  being  used.)  The  five  parameters  encompass 
eight  variables:  p,  C^,  a,  d  ,  u^,  g,  and  y  or  t.  The  results  of 
interest  are  shown  in  Figures  2-4  through  2-11  for  which  the  values  of 
Cjj  =  0.5  and  K  =  0.001  are  used.  The  value  of  =  0.5  is  a  good  approxi¬ 
mation  since  the  particle  Reynolds  number  is  generally  very  high,  and 
K  =  0.001  is  typical  for  blast  wave  conditions.  In  Figures  2-4  through 
2-11,  the  subscript  e  denotes  conditions  at  the  edge  of  the  boundary 
layer  (y^  =  4.66)  during  ascent,  m  at  the  maximum  height  reached  by  the 
particle,  and  f  at  the  point  where  the  particle  has  returned  to  the  edge 
of  the  boundary  layer  during  descent.  The  normalized  particle  velocity 
components  at  the  edge  of  the  boundary  layer,  which  are  obtained  by  assum¬ 
ing  drag  and  gravity  to  be  negligible  compared  to  lift  inside  the  boundary 
layer,  are  shown  in  Figure  2-4.  The  maximum  height  reached  by  the  particle 
is  shown  in  Figures  2-5  and  2-6  where  different  normalization  factors  are 

used.  From  Figure  2-5,  it  is  seen  that  the  particle  cannot  travel  far 

2 

above  the  boundary  layer  unless  Ku  is  large  (for  fixed  a  and  6).  It  can 

also  be  deduced  from  Figure  2-5  that  there  is  some  optimum  values  of  the 

ratio  of  particle  size  to  boundary  layer  thickness  for  aerodynamic  lift. 

2 

In  Figure  2-6,  the  nondimens ional  group  ^8  /ku^  is  employed  since 

it  can  be  shown  from  Equations  (2-26)  and  (2-31)  that  2g  (y^-y^) /Ku^^  ->■  1 
as  0  and  a/6  ->  0.  The  particle  vertical  velocity  at  y^  is  zero, 

while  the  horizontal  velocity  at  y^  is  shown  in  Figure  2-7.  The  velocity 

components  of  the  particle  when  it  has  returned  to  the  edge  of  the  boundary 
layer  are  shown  in  Figures  2-8  and  2-9  from  which  it  can  be  deduced  that 
the  particle  returns  to  the  surface  almost  horizontally.  The  nondimen- 
sional  total  horizontal  distance  traveled  by  the  particle  and  total  time 
of  particle  stay  in  the  gas  are  shown  in  Figures  2-10  and  2-11,  respec¬ 
tively  . 

It  was  indicated  earlier  that  particles  for  which  lift  is  not  sub¬ 
stantially  greater  than  drag  and  gravity  inside  the  boundary  layer  will 
not  achieve  a  substantial  vertical  velocity  at  the  edge  of  the  boundary 
layer.  It  is  interesting  to  observe  what  this  condition  implies. 
Considering  the  order  of  magnitude  of  the  various  terms  in  Equation  (2-10) 

I 
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1/2 

and  using  the  approximate  result  of  ^  K  u,  it  can  be  shown  that  the 

following  condition  must  hold  in  order  for  the  particle  to  trave]  far 
above  the  boundary  layer: 


1/2 

K  '  C  6 
1  >>  - ^ —  + 


_SlL 


(2-33) 


Ku 


1/2 

The  ratio  K  6/a,  defined  as  s  by  Equation  (2-24),  may  be  considered 

to  be  a  measure  of  the  relative  magnitudes  of  vertical  drag  and  lift. 

2  -1 

The  ratio  g6/Ku^  ,  which  is  equal  to  a  times  the  factor  r  given  by 
Equation  (2-28),  is  a  comparison  between  the  effects  of  gravity  and  lift. 
Furthermore,  the  condition  given  by  Equation  (2-15)  must  be  satisfied  in 
order  for  a  particle  to  be  lifted  off  a  surface. 

The  above  results  were,  obtained  by  consiaering  the  particle  to  be 
lifted  by  an  aerodynamic  force.  It  may  be  of  inr^rest  to  estimate  the 
velocity  of  a  particle  as  it  leaves  the  ground  if  the  particle  is  partly 
imbedded  in  the  soil  initially  where  the  soil  pressure  is  higher  than  the 
fluid  pressure.  For  this  purpose,  consider  the  soil  pressure  to  be  p^  and 
the  fluid  pressure  to  be  p.  Let  the  initial  distance  between  the  particle 
center  and  the  soil  surface  be  denoted  by  y^.  The  equation  for  the  particle 
motion  is  given  by 


dv  3(p 


P  dy 


-  P)  y2 

-  (1  „  jL-  ) 

4p  a  2 

s  a 


(2-34) 


For  the  boundary  condition  of  v^  =  0  at  y  =  y^,  the  solution  is 


V  = 

P 


2  P  a 
s 


y  -  Yi  - 


3  3^ 

y  "  y^  ' 

3a^  / 


1/2 


(2-35) 


The  maximum  v  is  obtained  at  y  =  a  and  for  y,  =  0  and  is  equal  to 

?/2  ^ 

[(p  -  p)/p  ]  .  If  (p  -  p)  is  taken  to  be  the  free-stream  dynamic 

®  ®  2  ^  1/2 
pressure,  pu  /2,  the  maximum  v  is  (p/2p  )  u  .  This  v  is  in  the 
e  1/2  P  s  e  p 

neighborhood  of  (3p/8p^)  u^  which  is  tlie  estimated  maximum  achievable 

due  directly  to  lift. 
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FIGURE  2-5  MAXIMUM  PARTICLE  HEIGHT 
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FIGURE  2-8  PARTICLE  VERTICLE  VELOCITY  AT  IMPACT 


■H 


CONFIDENTIAL 


16118-6002-R7-00 


3.  SAMPLE  TRAJECTORIES 

The  analysis  of  the  previous  section  shows  analytically  the  effects  of 
the  various  parameters  on  particle  trajectory  for  a  selected  lift  coefficient, 
A  uniform  steady  flowfleld  with  constant  boundary  layer  thickness  was  assumed. 
In  the  present  section,  the  particle  trajectories  in  a  blast  flowfleld  tran¬ 
sient  and  nonuniform)  which  includes  the  boundary  layer  properties  obtained 
in  Part  I  of  this  study  are  investigated.  The  equations  are  taken  as: 


dx 

dt  '^p 


V 

dt  p 


du 

^  =  F  (u-Up) 


dv 


«  F(v-v  )  +  K(u-u  )|^  -  g 
dt  p  P  oy 


dT 

-  G(T-I  > 


0 


da 
dt 

dT 

_ I 

dt 


c_„a 


for  T  <  T 
P  pv 


for  T  <  T 
P  PV 


for  T  »  T 
P  pv 


for  T  «  T 
P  pv 


(3-1) 

(3-2) 

(3-3) 

(3'4) 

(3-5) 

(3-6) 

(3-7) 

(3-8) 


where  K  •=  3p/8pg  and 


F  =  KVCjj/a 


G  =  Kk^TNu/pCp^a  T^ 


1/2 


(3-9) 

(3-10) 

(3-11) 
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In  the  above  equations,  F  and  G  are  factors  for  drag  and  convective  heat 
transfer,  respectively;  c  is  the  specific  heat  of  the  particle,  T  is  the 
vaporization  or  ablation  temperature,  L  is  the  heat  of  ablation,  and  is 
the  fluid  conductivity  at  the  ambient  temperature  T^.  The  term  involving 
3u/3y  in  Equation  (3-4)  is  used  to  account  for  lift.  The  effects  of  particles 
on  the  fluid  flowfield  are  neglected.  Since  mass  and  heat  transfers  are  only 
of  secondary  consideration  in  this  study,  the  effects  of  ablation  on  the 
drag  and  heat  transfer  coefficients  are  also  neglected.  Furthermore,  radia¬ 
tion  is  not  considered.  The  drag  coefficient  and  Nusselt  number  Nu  are 
approximated  by 


C 


D 


Ih. 

Re 


+  0.5 


(3-12) 


Nu  =  2  +  0.459  Re 


0.55 


(3-13) 


where  Re  is  the  Reynolds  lumber  based  on  particle  diameter. 

For  the  fluid  piuperties,  the  vertical  velocity  v  is  taken  as  zero.  The 
axial  velocity  u,  temperature  T,  and  density  p  are  approximated  by  u  =  u^ 
d-e"^'^^),  T  =  T^  +  (T^-T^)  (l-e“'’'^*'t)^  p  =  p^(T^/T),  respectively,  where 
the  subscript  e  denotes  properties  for  inviscid  flow  which  are  computed  using 
the  Taylor-Sedov  strong  shock  solution.  The  parameters  n,  I,  and  are 
taken  from  the  laminar  boundary  layer  integral-exponential  solution  given 
in  Part  I  of  this  report,  and  T^  is  the  wall  temperature.  The  particle  is 
considered  to  be  initially  at  rest  on  a  flat  surface. 

From  the  numerical  results,  it  is  found  that  a  particle  of  given  size 
is  not  subject  to  direct  lift  in  a  region  extremely  close  to  the  shock  front 
where  the  boundary  layer  is  thin  nor  in  the  region  fa;,  away  from  the  shock 
front  where  the  boundary  layer  is  thick.  This  is  in  agreement  with  the 
results  of  Section  2.  It  is  also  found  that  direct  lift  does  not  loft  a 
particle  to  sufficient  height  for  suspension  except  for  a  very  short  time 
after  explosion  when  the  free-stream  velocity  is  extremely  high.  However, 
small  particles  that  are  lifted  at  early  times  vaporize  soon  after  leaving  the 
boundary  layer  because  of  the  high  gas  temperature.  Most  of  the  particles 
that  get  lifted  and  do  not  ablate  are  confined  within  a  few  inches  above  the 
ground.  They  can  achieve  a  vertical  velocity  of  only  a  small  fraction  of 
the  free-stream  gas  velocity;  however,  they  quickly  attain  a  high  horizontal 
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velocity,  which  is  comparable  to  the  free-stream  velocity,  with  which  they 

return  to  the  ground.  This  substantiates  previous  assumptions  and  observa- 
12  3 

tions  ’  ’  that  particles  are  not  directly  lifted  into  suspension  by  a  free- 
stream,  Instead,  a  particle  may  bounce  off  the  ground  a  few  times  and  achieve 
a  higher  velocity  with  each  subsequent  bounce  until  the  velocity  is  suffi¬ 
ciently  high  to  loft  the  particle  or  bounce  off  another  particle  into  sus¬ 
pension,  However,  lift  does  aid  particle  entrainment  by  placing  the  particle 
in  the  higher  velocity  region  of  the  boundary  layer  to  facilitate  the  accel¬ 
eration  of  the  particle  by  the  gas  stream.  Also,  direct  lift  can  be  partially 
responsible  for  a  significant  amount  of  dust  having  bouncing  motions,  known 
as  saltation,  within  a  short  distance  from  the  ground.  Another  mechanism 
responsible  for  saltation  is  horizontal  drag  which  causes  a  particle  to  roll 
or  slide  along  a  surface  until  the  particle  accelerates  to  sufficient  speed  j 
to  leave  the  surface  upon  impacting  another  particle.  ' 

The  above  paragraph  describes  the  particle  trajectories  qualitatively. 

Some  numerical  results  are  presented  in  Figures  3-1  through  3-5  for  illus¬ 
tration.  The  blast  properties  correspond  to  a  1.0  megaton  spherical  surface 

3 

explosion.  The  physical  properties  used  are:  “  530'’R,  =  0.076U/ft  , 

C  =0.24  BTU/lb-®R,  h*,  =  1.23  X  10“^  Ib/ft-sec,  =  4.1  x  lO"^  BTU/sec-ft- 
'“R,  Y  =  1.4,  P  =  145  Ib/ft^,  c  =0.2  BTU/lb-^R,  T  =  6400  ®R,  and  L  = 

3700  BTU/lb.  Figure  3-1  shows  the  trajectories  (y  vs  x)  of  particles  of 
various  sizes;  the  particles  have  been  assumed  to  be  initially  at  rest  and 
get  lifted  at  0.01  second  after  explosion  and  at  a  distance  of  10  feet  behind 
the  shock.  Note  that  drastically  different  scales  have  been  used  for  the  x 
and  y  coordinates.  The  smaller  particles  become  ablated  by  the  hot  gas  soon 
after  they  leave  the  boundary  layer  and  are  indicated  in  Figure  3-1.  The 
larger  particles  are  lifted  to  higher  altitudes  (there  is  an  increase  in 
boundary  layer  thickness  as  the  particle  is  being  lifted),  but  there  is  a 
cutoff  diameter  (roughly  20,000  microns  for  the  conditions  of  Figure  3-1) 
above  which  a  particle  does  not  get  lifted  at  all  because  its  weight  becomes 
larger  than  the  lift  at  the  ground  position. 

Figures  3-2  to  3-4  are  similar  to  Figure  3-1  except  that  the  time  of  0,1 
second  after  explosion  is  used.  Also,  the  locations  at  which  the  particles 
become  lifted  are  1  foot,  10  feet,  and  100  feet  behind  the  shock  for  Fig¬ 
ures  3-2,  3-3,  and  3-4,  respectively.  At  0,1  second  after  explosion,  it  is 

*This  is  calculated  based  on  1.0  megaton  free  air  burst  which  corresponds 
roughly  to  1/2  megaton  surface  burst. 
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found  that  the  particles  do  not  reach  ablation  temperature  (except  for  small 

particles  at  locations  beyond  100  feet  from  the  shock) ,  The  initial  free- 

stream  gas  velocity  u^^  and  the  initial  velocity  boundary  layer  thickness 

y^jl^  (u  •  0.99  u^  at  y  *  y^)  »  along  with  the  velocity  components  v^^  and 

when  the  particles  impacts  the  ground  and  the  length  of  time  tt  in  which 

the  particle  stays  in  the  gas,  are  tabulated  in  the  figures.  In  Figure  3-2, 

the  location  of  the  shock,  x  >  when  the  particle  returns  to  the  ground  and 

sg 

the  velocity  boundary  layer  thickness,  indicated  by  y  ,  experienced  by  the 
particle  when  it  returns  to  the  ground  are  also  indicated.  It  can  be 
observed  that  a  particle  accelerates  more  rapidly  towards  the  ground  as  it 
returns  to  the  boundary  layer  region;  this  is  due  to  a  negative  lift  which 
develops  since  u^  becomes  greater  than  u  inside  the  boundary  layer.  How¬ 


ever,  V  is  still  very  small  compared  to  u  , 
Pg  Pg 

strike  the  ground  almost  horizontally. 


and  a  particle  is  expected  to 


From  Figures  3-1  to  3-A,  it  is  seen  that  the  particles  are  lifted  not 
more  than  a  few  inches  above  the  ground  except  for  early  time  after  explo¬ 
sion  or  from  regions  of  large  boundary  layer  thickness.  It  is  also  seen 
from  the  figures  that  for  a  given  boundary  layer  thickness,  there  is  some 
optimum  particle  size  for  entrainment.  Particles  can  be  lifted  to  higher 
altitudes  if  their  diameter  is  roughly  equal  to  the  boundary  layer  thickness, 
i.e.,  D/y^^  ss  1.  It  is  interesting  to  note  that  the  analysis  of  Section  2 
predicts  a  similar  result:  in  Figure  2-5,  y^  is  largest  for  a/6  =1.5;  since 
a  =  D/2  and  6  =  y^/4.6,  the  maximum  y^  corresponds  to  D/y^  =  0.7. 

After  returning  to  the  ground  at  high  speed,  the  particles  are  expected 
to  rebound  from  the  ground.  The  rebound  speed  and  direction  depend  on  the 
soil  properties  such  as  elasticity  and  surface  smoothness.  Some  rebound 
trajectories  are  illustrated  in  Figure  3-5.  The  position  from  which  a 
particle  rebounds  corresponds  to  the  end  of  the  particle’s  lift  trajectory 
shown  in  Figure  3-3.  The  rebound  speed  is  taken  to  be  the  impact  speed, 
and  the  rebound  angle  is  taken  to  be  0,  45,  or  -45  degree  from  the  normal 
direction  to  surface.  As  expected,  the  larger  particles  can  rebound  to  much 
greater  heights  than  the  smaller  particles  because  of  lower  drag  forces. 
Figure  3-5  also  shows  that  the  particle  trajectories  are  similar  for  all 
three  rebound  angles  chosen  for  illustration.  The  end  of  each  trajectory 
corresponds  to  2.0  seconds  after  explosion.  Thus,  a  particle  that  is  lifted 
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off  a  purfacc  may  travel  only  a  short  vertical  distance  and  return  quickly 
to  the  ground;  but  after  first  bounce,  it  may  remain  in  the  gas  for  several 
seconds  which  is  sufficiently  long  for  the  particle  tc  be  carried  aloft  by 
the  rising  fireball  thermal,  Furthermore,  the  larger  particles  can 
achieve  a  height  of  several  hundred  feet  simply  by  rebound. 
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FIGURE  3-2  PARTICLE  TRAJECTORY  FOR  A  1.0  MEGATON  EXPLOSION  WITH  INITIAL  TIME  0.1  SECOND 
AFTER  EXPLOSION  AND  INITIAL  DISTANCE  1  FOOT  BEHIND  SHOCK 


16118-6002-R7-00 


CONFIDENTIAL 


FIGURE  3-4  PARTICLE  TRAJECTORY  FOR  A  I.O  MEGATON  EXPLOSION  WITH  INITIAL  TIME  0.1  SECOND 
AFTER  EXPLOSION  AND  INITIAL  DISTANCE  100  FEET  BEHIND  SHOCK 


FIGURE  3-5  PARTICLE  REBOUND  TRAJECTORY 
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4.  EROSION  ESTIMATES 

In  order  to  analyze  rigorously  the  amount  of  soil  erosion  due  to  aero¬ 
dynamic  entrainment,  it  is  necessary  to  solve  the  coupled  two-phase  (particles- 
gas)  Navier-Stokes  equations.  The  development  of  such  a  solution  for  the 
blast  flowfield  appears  prohibitive  at  present,  however,  since  even  the 
simple  case  of  erosion  over  a  flat  plate  has  not  been  solved.  Here,  only 
an  extremely  simple  model  is  used  to  estimate  erosion.  The  model  is  similar 
to  the  one  first  employed  by  Roberts^  in  regard  to  soil  erosion  by  a  lunar 
landing  vehicle,  and  it  can  be  written  as 

£  dm  * 

fu  -Tr  =  -  T  (4-1) 

at  o 

where  f  is  interpreted  as  a  dimensionless  average  momentum  of  the  entrained 
particles,  u  is  a  gas  velocity,  dm/dt  is  the  rate  of  mass  entrained  per  unit 
area,  t  is  the  shear  stress  on  the  surface  without  dust  entrainment,  and 
T  is  the  reduced  shear  stress  due  to  the  presence  of  dust  entrainment. 

Equation  (4-1)  states  that  the  impulse  due  to  shear  [(’’^q  -T*)dt]  is  con¬ 
verted  into  particle  momentum  with  the  particle  velocity  being  fu.  The 
value  T*  represents  the  soil's  resistive  shear  stress  to  erosion,  and 
depends  on  the  soil  properties  such  as  particle  size,  cohesiveness, 
etc.  The  value  of  f  is  difficult  to  analyze;  in  fact  it  partly  depends  on 
the  definition  of  mass  erosion.  If  only  those  particles  that  reach  high 
altitude  are  considered,  then  f  should  be  close  to  1.0  since  only  the  parti¬ 
cles  with  initial  velocities  that  are  comparable  to  the  gas  velocity  can  be 
lifted  to  significant  height.  On  the  other  hand,  a  great  many  particles 
that  are  lifted  off  a  surface  according  the  results  of  Section  2  are  confined 
within  a  short  distance  from  the  surface;  these  particles  have  lower  veloc¬ 
ities  as  they  leave  the  boundary  layer  and  f  should  be  substantially  less 
than  1.0.  In  this  section,  f  is  taken  as  an  arbitrary  constant  and  t*  is 
neglected  when  it  is  compared  to  t^.  The  total  mass  M  eroded  is  then  given 
by 


M  .  7  /  /  —  dAdt  (4-2) 

Using  the  results  of  Part  I  of  this  study  (laminar  boundary  layer  integral 
polynomial  solution)  for  and  Sedov's  strong  shock  expression  for  u,  it  is 
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found  that  for  a  spherical  explosion 


M  =  /l±i 

f  \  2  p„/ 


(4-3) 


where  y  is  the  ratio  of  the  gas  specific  heats,  C  is  a  constant  which  depends 

on  the  explosion  strength,  and  u  ,  p^,  and  denote  the  gas  viscosity, 

pressure,  and  density,  respectively,  at  ambient  conditions.  Since  C  varies 

with  the  explosion  energy  E  to  the  0.2  power.  Equation  (4-3)  shows  M  to  vary 

with  E  to  the  0.6  power.  Also,  M  varies  with  time  t  to  the  0.7  power.  Here, 

t  is  bounded  by  the  time  when  the  shock  is  no  longer  strong  or  when  the  soil 

1  exceeds  t  .  Consider  a  1.0  megaton  explosion  (i.e.,  C  =  3800  ft/sec  '  ) 

in  ambient  air  and  take  t  =  1.0  second.  Then  M=  35  f  tons.  Hence  if  f  is 

close  to  1.0,  M  is  very  small  (compared  to  mass  ejected  due  to  cratering). 

However,  if  the  particle  velocity  is  taken  as  the  vertical  velocity  due  to 

1/2 

boundary  layer  lift,  then  f  a  K  '  where  K  =  0.001.  This  yields  M  =  1,100 
tons  for  laminar  boundary  layer  shear. 

It  may  be  mentioned  that  the  total  impulse  delivered  to  the  surface 
due  to  aerodynamic  shear  is  found  to  be  (for  laminar  boundary  layer) : 


J  /  T  dAdt  =  1.82  — ) 

t  A  o  \  2  p^7 


Poo^  t 


(4-4) 


For  a  1.0  megaton  explosion  in  ambient  air  and  t  =  1.0  second,  the  total 

g 

impulse  is  equal  to  4.15  x  10  Ib-ft/sec.  If  the  mass  eroded,  M,  is  taken 
to  be  equal  to  this  impulse  divided  by  some  effective  average  particle  veloc¬ 
ity,  then  M  is  proportional  to  E*^'®  (since  C  ~  E*^’^)  and  t^’^;  whereas 
in  Equation  (4-3),  M  is  proportional  to  E*^’^  and  t^'^. 

The  numbers  given  in  the  two  preceeding  paragraphs  are  for  laminar  flow. 
For  turbulent  flow,  the  values  for  M  is  expected  to  Increase  by  at  least  one 
order  of  magnitude.  Although  the  assessment  of  soil  erosion  behind  the 
turbulent  blast  wave  should  probably  be  performed  only  after  the  turbulent 
boundary  layer  solution  is  obtained,  it  may  be  interesting  to  make  a  predic¬ 
tion  by  employing  a  simple  model  and  available  experimental  data.  Consider 
the  model 


(4-5) 
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where  f!  is  a  constant.  Equation  (4-5)  states  that  the  mass  erosion  rate 
is  proportional  to  the  free-stream  velocity.  This  model  is  similar  to  that 
employed  in  Reference  6,  except  particle  size  and  terminal  particle  velocity 
are  not  considered  and  a  factor  P /p^  is  inserted  to  account  for  the  varying 
density  behind  a  blast  wave.  The  mass  eroded,  M,  then  becomes 

M  =  6  Pg  ^  j  dAdt  (4-6) 

A  t 

Using  the  Sedov  solution  for  p  and  u  for  a  spherical  explosion,  it  is 
found  that 


M  =  0.55BP  (4-7) 

s 

3  0  6  1  2 

Here,  M  is  proportional  to  C  (or  E  )  and  t  *  .  Refering  to  the  experi¬ 
mental  data  of  Reference  6,  one  may  suspect  that  B  is  far  from  being  a  con- 

—5  -3 

stant  and  may  vary  between,  oay,  10  and  10  .  In  fact,  6  for  loose  soil 

is  about  one  order  in  magnitude  higher  than  that  for  compact  soil.  For  a 

given  soil,  B  varies  with  u.  One  may  suspect  that  the  scatter  in  data  for 

B  is  at  least  partly  due  to  the  data  correlation  in  which  the  direct  effects 

of  boundary  layer  growth  are  not  taken  into  account.  Taking  B  =  6  x  10  ^ 

1  5 

and  Pg  =  145  Ib/ft  ,  it  is  found  that  M  =  1.32  x  10  tons  for  a  1.0  megaton 
explosion  and  t  =  1,0  second. 


(Reverse  of  Page  is  Blank) 
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5.  DISCUSSION  AND  CONCLUSIONS 

Using  a  lift  coefficient  of  C  =  a(u  -  u  )  ^  (du/dy)  for  boundary 

L  p 

layer  aerodynamic  lift,  it  is  found  in  the  present  study  that  a  particle 

1/2 

initially  at  rest  can  attain  a  vertical  velocity  of  only  K  times  the 

free-stream  gas  velocity  where  K  is  equal  to  3p/8p  .  A  similar  vertical 

s 

velocity  is  attained  if  a  particle,  is  initially  Imbedded  in  the  soil  and 
if  the  pressure  difference  between  the  soil  and  the  free-stream  is  equal  to 
the  dynamic  pressure.  Also,  spinning  of  a  particle  may  generate  a  Bernoulli 
force  for  which  the  lift  coefficient  is  similar  to  that  due  to  shear. 

The  analytical  results  of  Section  2,  which  are  for  constant  flow  condi¬ 
tions,  may  be  applied  to  obtain  rough  estimates  of  particle  motion  due  to 
boundary  layer  lift  behind  blast  waves.  The  reason  for  this  is  that  a 
particle  stays  in  the  gas,  especially  inside  the  boundary  layer,  for  only 
a  short  time  within  which  the  fluid  properties.  Including  the  boundary 
layer  thickness,  experienced  by  the  particle  do  not  change  drastically. 

For  given  flow  conditions  and  a  given  particle  size,  the  aerodynamic 
lift  is  not  sufficient  to  overcome  the  particle  weight  to  lift  the  particle 
if  the  boundary  layer  thickness  is  either  too  small  or  too  large.  Particles 
that  are  more  readily  entrained  are  those  with  diameters  of  the  order  of 
the  boundary  layer  thickness. 

In  general,  particles  are  not  lifted  more  than  a  few  inches  above  the 
ground  except  for  early  times  after  explosion  when  the  gas  velocity  is 
high  and/or  at  locations  far  away  from  the  shock  (but  not  too  far)  where 
the  boundary  layer  is  thick  and  the  gas  velocity  and  density  are  still 
high.  However,  at  ear.ly  times  after  explosion  or  at  distances  near  the 
explosion  point,  the  temperature  is  extremely  high  and  the  lifted  particles 
(especially  the  small  ones)  soon  vaporize. 

As  a  particle  is  lifted  off  a  surface,  it  attains  a  very  high  hori¬ 
zontal  velocity.  Although  it  reaches  a  vertical  height  of  only  a  few 
inches  and  stays  in  the  gas  for  only  a  fraction  of  a  second,  it  may  travel 
a  horizontal  distance  of  a  few  hundred  feet  before  returning  to  the  ground. 
It  impacts  the  ground  about  horizontally  at  a  high  speed.  The  height 
reached  by  the  particle  after  rebounding  from  the  ground  is  still  small  for 
small  particles  (say  less  than  100  microns  in  diameter)  because  of  drag 
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effects,  but  the  height  can  be  on  the  order  of  a  few  hundred  feet  for 
large  particles  (say  lOOO's  of  microns  in  diameter).  The  particles  after 
rebounding  can  remain  in  the  gas  for  several  seconds.  Thus,  although 
lift  may  not  directly  contribute  significantly  to  particle  density  in  a 
nuclear  cloud,  it  can  aid  the  particle  in  acquiring  a  high  rebounding 
velocity.  The  particles  after  rebound  can  stay  in  the  gas  sufficiently 
long  to  be  carried  to  great  heights  by  a  rising  nuclear  cloud. 

In  the  present  study  on  particle  lift,  the  vertical  component  of  the 
gas  velocity  has  been  neglected.  This  vertical  velocity  may  substantially 
increase  the  particle  vertical  velocity,  especially  for  small  particles,  ! 

I 

in  the  regions  of  thick  boundary  layers.  Turbulent  diffusion  also  aids 
entrainment  of  small  particles. 

I 

The  amount  of  soil  erosion  due  to  laminar  boundary  layer  entrainment 
is  estimated  to  be  small,  say  about  one  thousand  tons  for  a  one  megaton  j 

explosion.  However,  for  turbulent  boundary  layers,  an  increase  by  two  i 

I 

orders  of  magnitude  appears  possible.  I 


] 
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